Discovery of transient infrared emission from dust
heated by stellar tidal disruption flares
Stars that pass within the Roche radius of a supermassive black hole will be tidally disrupted, yielding a sudden injection of gas close to the black hole horizon which produces an electromagnetic flare. A few dozen of these flares have been discovered in recent years, but current observations provide poor constraints on the bolometric luminosity and total accreted mass of these events. Using images from the Wide-field Infrared Survey Explorer (WISE), we have discovered transient 3.4 m emission from several previously known tidal disruption flares. The observations can be explained by dust heated to its sublimation temperature due to the intense radiation of the tidal flare. From the break in the infrared light curve we infer that this hot dust is located pc from the supermassive black hole. Since the dust has been heated by absorbing UV and (potentially) soft X-ray photons of the flare, the reprocessing light curve yields an estimate of the bolometric flare luminosity. For the flare PTF-09ge, we infer that the most likely value of the luminosity integrated over frequencies at which dust can absorb photons is erg s, with a factor of 3 uncertainty due to the unknown temperature of the dust. This bolometric luminosity is a factor larger than the observed black body luminosity. Our work is the first to probe dust in the nuclei of non-active galaxies on sub-parsec scales. The observed infrared luminosity implies a covering factor % for the nuclear dust in the host galaxies.
Studying the inner parsec of galaxies with imaging observations is notoriously difficult due to the very high concentration of stars and the finite resolution of our telescopes. Fortunately, reverberation of time-variable signals can provide information on scales much smaller than can be spatially resolved. The response of broad emission lines to variability in the continuum light of active galactic nuclei (AGN), for example, yields a scale for their emission region of light days to light years (Kaspi et al., 2000). Hot dust in a torus around the AGN disk can reprocess the optical/UV flux from the disk, re-emitting it at a typical temperature K (Barvainis, 1987). Indeed the infrared (IR) light curve is observed to lag the optical light curve and implies a torus size of pc (with the hard X-ray luminosity of the AGN, Koshida et al. 2014).
Similar to reverberation mapping of AGN, gas at the centers of inactive galaxies can be studied via light echoes produced after the central supermassive black hole tidally disrupts a star. A few dozen tidal disruption flares (TDFs) candidates have been discovered in the last decade, first using X-ray (e.g., Komossa & Bade, 1999; Esquej et al., 2008; Saxton et al., 2012) and UV satellites (Gezari et al., 2006, 2009), and more recently using optical surveys (e.g., van Velzen et al., 2011; Gezari et al., 2012; Chornock et al., 2014; Arcavi et al., 2014; Holoien et al., 2014, 2015).
Light echoes from TDFs could be detected by monitoring (broad) emission lines in optical spectra of these flares (Komossa et al., 2008). However the origin of observed TDF emission lines is unclear: the lines may not be due to photoionization of the pristine gas around the black hole, but may instead originate from within the photosphere of the stellar debris (Roth et al., 2015).
When a stellar tidal disruption results in the launch of a powerful relativistic jet (Zauderer et al., 2011; Bloom et al., 2011; Burrows et al., 2011; Levan et al., 2011; Cenko et al., 2012), the radio light curve of this jet can be used to estimate the density of the circumnuclear gas as a function of distance to the black hole (Berger et al., 2012). Yet this inference requires the assumption that the microphysical parameters of the jet (e.g., the fraction of energy in magnetic fields) remain unchanged as it expands.
A new, and hitherto untested, use of TDFs to study the nuclei of inactive galaxies is transient IR emission from dust heated by the intense radiation field of the flare (Lu et al., 2016).
The typical scales relevant to IR reprocessing can be estimated from existing data and dust grain physics. The few dozen TDF candidates discovered so far are bright () month-long flares that can be described by a single-temperature black body at optical/UV frequencies, K (van Velzen et al., 2011; Gezari et al., 2012), or soft X-ray frequencies, keV (Komossa & Bade, 1999; Miller et al., 2015). This high luminosity implies any dust near the black hole will quickly evaporate. The dust temperature is set by the flux absorbed by the grains and thus decreases with increasing distance to the black hole. For graphite grains with a radius of 0.1 m, we find
(Barvainis, 1987; Waxman & Draine, 2000). Here is the dust temperature and is the flare luminosity integrated over the frequencies where dust absorbs. We normalized to the expected sublimation temperature of the dust (discussed in detail in Sec. 2). If 1% of the radiated TDF energy ( erg) is reprocessed and re-emitted over one year, the expected IR luminosity is or % of the typical galaxy luminosity in the -band.
To search for a dust reprocessing signal from TDFs, we used Wide-field Infrared Survey Explorer (WISE; Wright et al., 2010) images at 3.4 and 4.6 m. Besides a very stable photometric performance, an advantage of using WISE observations for this search is a factor of 2 reduced host galaxy flux at 3.4 m compared to the -band.
Our work is the first111 Two preprints on TDF dust reprocessing (Dou et al., 2016; Jiang et al., 2016) appeared a few days after the preprint of this manuscript was posted to the arXiv. Our discovery of dust reprocessing in WISE data was made public in van Velzen et al. (2015). to present transient IR emission from thermal (i.e., non-relativistic) TDFs, showing that dust heated by a TDF is detectable. These observations thus open a new wavelength regime for the study of TDFs and provide a new probe to study dust in galaxy centers.
Prior to the work presented here, IR observations of TDFs have been focused mainly on the relativistic TDF Swift J1644+57 (Bloom et al., 2011; Zauderer et al., 2011; Burrows et al., 2011; Levan et al., 2011), whose bolometric output is dominated by synchrotron emission from a relativistic jet. The -band and Spitzer 3.6 and 4.5 m emission of Swift J1644+57, observed about 300 days after discovery (Levan et al., 2016; Yoon et al., 2015), appear to fall on the Rayleigh-Jeans tail of a blackbody with a temperature K (Lu et al., 2016). This thermal component could be due to a stellar explosion or an accretion disk (Levan et al., 2016), but the late-time IR emission of Swift J1644+57 is not consistent with dust reprocessing.
The TDF candidate SDSS J0952+2143, identified in the Sloan Digital Sky Survey (SDSS; York et al., 2000) spectroscopic sample based on its extreme coronal emission lines (Komossa et al., 2008), has been followed-up with Spitzer mid-IR spectroscopic observations to search for signs of dust heated by the flare. These observations were obtained about 4 years after the peak of the optical light curve (Palaversa et al., 2016). WISE observations obtained 6 years after the peak suggest this Spitzer flux has faded by a factor of two (Dou et al., 2016), providing strong evidence for transient IR emission due to dust reprocessing.
The outline of this paper is as follows. We first present the necessary theoretical background for dust reprocessing in the context of TDFs in Section 2. We then discuss the steps of our data analysis in Section 3: sample selection (Sec. 3.1), co-addition of WISE images (Sec. 3.2), and photometry (Sec. 3.3). In Section 4, we explain how the parameters of our dust reprocessing model are estimated from the IR light curves. We discuss the results in Section 5 and close with a list of our conclusions. We adopt a flat CDM cosmology with and km s
2. Dust reprocessing for tidal disruption flares
Because hydrogen ionization can remove photons before they heat the dust, we first estimate the magnitude of this effect and find that both H photoionization and H photodissociation absorb only a very small fraction of the TDF light.
In the initial stages of the TDF, its optical/UV continuum rises steeply toward a maximum. As these first photons travel outward through surrounding gas, they ionize its H atoms. The spectra of most TDF are reasonably well-fit by blackbodies of constant temperature K (Gezari et al., 2009; van Velzen et al., 2011; Arcavi et al., 2014; Chornock et al., 2014; Holoien et al., 2015). In such a case, the ratio of the number of ionizing photons radiated over the course of the flare to the number of nearby H atoms is
where is the integrated optical/UV energy in units of erg, is the distance of the ionized region in units of 0.1 pc, and is the volume density of H atoms in units of cm. Thus, the flare is more than capable of ionizing all the atoms out to nearly a parsec unless the density is cm. The availability of ionizing photons diminishes for lower blackbody temperatures, but drops by only a factor of 30 if the temperature is as low as K. If the gas is initially molecular rather than atomic, the relevant photons are in the 7.5–13.6 eV range, which are even more abundant in the light radiated by the flare than those above the ionization edge, so H molecules will be broken into separate atoms even more quickly than individual H atoms are ionized. Moreover, because the recombination time is year, these atoms stay ionized for at least several years unless cm.
The flare is similarly able to evaporate the nearer dust grains. The binding energy per atom is 5.7 eV for astronomical silicate or 6.8 eV for graphite (Guhathakurta & Draine, 1989). If the total number density of C, N, O, and Si atoms is (Asplund et al., 2009) and half of them are in grains, the ratio of radiated flare energy to grain binding energy is
This order of magnitude estimate shows that sublimating all dust grains within 0.1 pc would remove only a tiny fraction of the TDF photons.
While there is enough energy in the optical/UV flare to vaporize the grains out to parsec distances, the time required for grains to sublimate can be significant. Following Guhathakurta & Draine (1989), we find the grain survival time,
where is the grain temperature, the factor 42.747 applies to graphite grains, and we have normalized the grain size to 0.1 m.
As the dust evaporates, the IR flux fades rapidly. For a single dust grain,
with the emission efficiency. At the wavelength of our observations and for grains with m, (cf. Figure 4b in Draine & Lee, 1984). Since , we have
and we find the effective sublimation temperature () by equating one fourth of the sublimation time for a given grain (Eq. 4) to the duration of the TDF. The full width at half maximum (FWHM) of the optical TDF light curve of PTF-09ge is months, which yields K. At this temperature, silicate dust sublimates in hr, hence this grain type can survive only at larger radii from the black hole compared to graphite. As we will show below, a large radius leads to a lower IR luminosity because the reprocessed energy is emitted over a longer time. We therefore model the reprocessing signal using only graphite grains. The sublimation temperature is weakly dependent on grain size, e.g., for m, the sublimation temperature is 100 K lower.
Following estimates of the dust heated by AGN or GRBs (Barvainis, 1987; Waxman & Draine, 2000), we compute the radius of a dust shell as a function of its temperature by equating the rate at which the grains absorb heat,
to the rate at which dust radiates (Eq. 5):
Here is the absorption efficiency of the grains and is the temperature-averaged emission efficiency. The latter depends on the shape of the dust spectral energy distribution (SED). Using a modified black body spectrum, , with at the wavelength of our observations (Draine & Lee, 1984), the emission efficiency is
(Draine & Lee, 1984, Eq. 6.1) and we thus find the radius of the dust shell as function of the dust temperature, size and luminosity of the heat source:
Here and K. For , Eq. 10 yields the sublimation radius. Outside the sublimation radius, the grain temperature stays in close equilibrium with the incident radiation for the grain cooling time is extremely small:
As will be discussed below in Section 4, the shape of the IR light curve can be used to estimate the radius from the black hole where this emission originates. We therefore rewrite Eq. 10 to find the luminosity of the flare as a function of this radius:
From this expression we see that the flare luminosity inferred from the reprocessing light curve is sensitive to both the dust temperature and the dust grain size. The dust temperature can be measured using multi-frequency IR follow-up observations of TDFs; an accurate measurement of this temperature coupled with the flare duration would place a lower bound on the characteristic grain size.
For a typical distribution of grain sizes (; Weingartner & Draine 2001), the largest graphite grains determine the effective sublimation radius because they dominate the luminosity at 3 m:
The size normalization used in the expressions for and (i.e., m), is motivated by dust models for the extinction curve of the Milky Way and the Magellanic Clouds (Weingartner & Draine, 2001), which yield a power-law dust-size distribution with an exponential cutoff at .
The grains heated by the flare, but not evaporated, reradiate the energy absorbed in the near-IR. However, because they are distant from the black hole, there is a significant delay before the IR reaches a distant observer. To be specific, if the direction from the black hole to the observer is the polar axis of a system of spherical coordinates, the delay is given by
where is the radial coordinate of a particular dust grain and is its polar angle. In response to an isotropically-radiated optical/UV flare with light curve , surrounding material produces an IR light curve
where is the response function. For a spherical shell of material at radius responding linearly to the irradiating flux,
where is the marginal IR emissivity (when the dust radiates in the IR exactly the same energy it absorbs from optical/UV/X-ray flare light, the marginal IR emissivity is identical to the covering factor). The final value is obtained only for , otherwise the root of the function’s argument lies outside the range of . In other words, spherical shells of linear reprocessors generically produce square-wave response functions extending from to .
Our fiducial model is a thin and spherically symmetric shell of dust that is optically thin to its own IR emission. If , it predicts a square wave, independent of the detailed shape of the optical/UV light curve. In Section 4 we show that applying this transfer function to the TDF light curve provides a good description of the IR observations.
In the previous section we found that IR emission from dust heated by a TDF yields a signal that can be detected up to year after the peak of the flare. Below we present the details of our search for this signal using WISE images at 3.4 and 4.6 m.
We typically have two or three observations separated by 6 months in the period 2010–2012, and, thanks to the reactivation of WISE (Mainzer et al., 2014), four more images from 2014 through 2015. Our parent sample (discussed in Sec. 3.1) consists of optical TDFs. The flux of these targets was extracted from co-adds of the WISE exposures (discussed in Sec. 3.2) using a forced photometry method (discussed in Sec. 3.3). The uncertainty on this flux was computed by repeating our photometry on a set of reference sources with a flux similar to the target, and thus includes both statistical and systematic sources of uncertainty.
Note. – The third column, , lists the time of maximum light of the optical/UV TDF light curve (except for TDE2, which was observed post-peak and we use the first observation of this source). and denote the TDF black body temperature and luminosity, respectively. The last column, , lists the total -band magnitude of the host galaxy.
3.1. Parent sample of TDFs
To be able to use the WISE observations of 2014–2015 as a baseline, we restricted our sample to TDFs that occurred between 2007 and 2010, leaving five optical/UV-selected flares (X-ray TDFs are excluded from the sample because their light curves are so sparsely sampled that the time of maximum light is poorly constrained). The properties of our five targets are listed in Table 1, below we provide a brief review of these sources.
Three of the five TDFs in our parent sample were found in the Palomar Transient Factory (PTF; Law et al., 2009) data by Arcavi et al. (2014): PTF-09ge, PTF-09djl, and PTF-09axc. Arcavi et al. selected relatively bright optical transients, . Optical spectra of the flares were obtained within a few months of the peak of the flare, covering the wavelength range 300-1000 nm (in the observer-frame). These spectra were used to estimate the black body temperature of these events. The optical spectra imply a post-starburst (or E+A) classification for the host galaxies of these three flares (Arcavi et al., 2014).
One source in our sample, TDE2, was found in SDSS imaging observations by van Velzen et al. (2011). This flare was also detected in GALEX imaging, about 1 year after the optical peak. The black body temperature was measured using SDSS (, , , ) photometry.
One source in our sample, D23H-1, was found by Gezari et al. (2009) using the Galaxy Evolution Explorer (GALEX; Martin et al., 2005) Time Domain Survey (Gezari et al., 2013). This TDF was observed with Chandra, both 3 and 116 days after its peak; no source was detected in these X-ray follow-up observations, implying . The black body temperature of this flare was computed from -band, NUV and FUV photometry.
The optical spectrum of the host galaxy of D23H-1 revealed narrow H emission and a Balmer decrement. This implies ongoing star formation () and a modest amount of extinction in the host galaxy: (Gezari et al., 2009). Correcting for this extinction, Gezari et al. found a factor 10 increase to the TDF black body luminosity relative to the luminosity derived from the uncorrected SED. It should be noted that the Balmer decrement is a galaxy-averaged property and provides a measurement of dust extinction toward the starforming regions where the Balmer lines are produced. The light from a TDF, on the other hand, samples dust in a narrow pencil beam to the galaxy center, hence the dust extinction derived from the Balmer decrement may not apply directly to TDFs. For this reason, and because the host galaxy spectra of the other TDFs in our sample do not allow an extinction measurement via the Balmer decrement, we will use the uncorrected black body luminosity of D23H-1.
3.2. Co-addition of WISE images
The WISE satellite completed full sweeps of the sky in six months following great circles with the Sun at the center. There are about 15 orbits per day and most sources are observed 12 times in each scan (Wright et al., 2010). For the sources in our parent sample, we first investigated the WISE multi-epoch photometry catalog, which provides a profile flux measurement for each exposure in which a source is detected, finding significant variability for most of our targets. To obtain a higher signal-to-noise measurement, we co-added the individual exposures using the icore software package (Masci, 2013) — which is similar to the software used by the WISE team to create their “Atlas Images” (Masci & Fowler, 2009).
We produced co-adds of the individual “level 1b” images of the W1-band (3.4 m) and W2-band (4.6 m) — the longer-wavelength bands (W3 and W4) are not available after September 2010, when the mission ran out of cryogenic coolant (Mainzer et al., 2011). We excluded frames flagged for having low image quality. We adopted the default parameters of icore for WISE data and produced area-weighted co-adds with a pixel scale of 1.0” pix (compared to 2.75” pix for the input images). No “drizzling” was applied, i.e., the pixel size of the input images was kept fixed when mapped to the co-add frame.
To have a sufficient number of reference sources in the field (see next section) we used a co-add size of (about a factor 2 smaller than the WISE field of view). The depth of coverage (i.e., number of individual images contributing to a co-add pixel) over this area is constant at the 5–10% level; for a box of 10 pixels, the depth of coverage is constant at the 3–5% level.
For each co-add, we estimated the point-spread-function (PSF) by fitting a Gaussian ellipse to isolated stars (Jones et al., 2015); we typically use 5–10 stars for this fit. The residuals of the fit are added to the Gaussian profile to find a good approximation of the true PSF. The typical PSF FWHM is 6”. To enable relative photometry, we also co-added all individual exposures (i.e., all images from 2010 to 2015). Cutouts of these co-adds are shown in Fig. 1.
|name||# of refs.|
Note. – For each TDF in our parent sample, we list the number of reference sources that were used by our relative photometry method (second column) and the mean W1 magnitude (third column). The mean of the rms scatter of the W1 light curve of the reference sources, , measures the typical accuracy of our photometric method. The over the degrees of freedom (dof) is computed for the W1 light curve of each TDF, using a model that has a constant flux as a function of time. The last column lists the probability to find the observed .
3.3. Photometry on WISE images
After obtaining a set of co-adds centered on each TDF, the next step is to extract the flux from these images. Given the depth and resolution of our co-adds, we anticipate blending of sources needs to be accounted for to obtain accurate photometry. Fortunately, all our targets are in the SDSS footprint, hence we can use information from the higher-resolution SDSS imaging to guide our WISE photometry. Our approach is motivated by the success of the unWISE project (Lang, 2014; Lang et al., 2014), which used the measured SDSS source positions and star/galaxy profiles to construct a catalog of WISE photometry for 400 million SDSS sources.
To measure the flux of a given target, we place it at the center of a scene. We then measure the flux of all known SDSS sources in this scene by minimizing the difference between a model of these sources and the observed image using galfit (Peng et al., 2002). For sources classified as stars by the SDSS pipeline (Stoughton et al., 2002), this model is simply given by the PSF of the co-add. For SDSS galaxies, we use a deVaucouleurs profile with parameters measured by the SDSS pipeline. Except for sources in the scene that are 2 mag brighter than the main target, we keep the centroid of the model fixed at the location measured by SDSS. By simultaneously fitting for the amplitude of all objects in the scene, we properly account for any contributions from nearby sources to the flux of the target. To obtain a stable solution, we fixed the amplitude of sources that are 2 mag fainter than the target at the magnitude reported by unWISE or sources with a distance from the center of the scene that is larger than 3 times the PSF FWHM.
Our method of image co-addition leads to spatially correlated noise, hence the uncertainty obtained from fluctuations of the “sky” underestimates the true statistical uncertainty. To measure the uncertainty, we use a set of reference sources (both stars and galaxies) that are assumed to have a constant flux with time. We measure the flux of the reference sources using the galfit scene photometry described above. By comparing the observed reference source flux in each epoch, , to the flux in the co-add of all exposures, , we can estimate the uncertainty on the flux in each epoch by fitting a Gaussian distribution to . We also use the mean of this distribution to correct for any offsets in the zero point of the co-adds with respect to the zero point of the co-add of all exposures. The relative photometry light curves, as shown in Fig. 2, are thus given by
With running over the reference sources. The relative zero point offsets, , are found to be between 0.001 mag to 0.015 mag. To recover the absolute flux scale, we tied our measurement of the flux of the reference sources to the flux listed in the unWISE catalog.
Reference sources are selected in a range of mag from the WISE flux of the TDF, with the exception of the brightest target (PTF-09ge), for which we use mag to ensure sufficient reference sources are available. Reference sources that are not consistent with a constant flux are removed using a simple cut of . The final number of reference sources for each TDF is listed in Table 2.
Since each reference source is expected to have constant flux with time, the observed variability of the reference sources provides a measurement of the typical accuracy of the photometry for each target. We used Eq. 17 to find the relative photometry light curve of each reference source and computed the root-mean-square variability of this light curve, . We list the mean of in the third column of Table 2. For PTF-09ge, the brightest target in our sample, the photometric accuracy is 0.016 mag.
Three of the five TDFs in our sample show significant variability (i.e., the probability of a constant flux is less than 1%, see Table 2) in the W1 band: PTF-09ge, PTF-09axc, and D23H-1. The signal-to-noise of the WISE observations of the other TDFs in our sample is not sufficient to detect variability at the level observed for these three flares.
To find the difference flux, we use co-adds of 2014-15 to estimate the baseline host galaxy flux. A change in flux relative to this baseline is simply given by
The uncertainty on follows from adding the uncertainty of the baseline and the pre-2014 measurements in quadrature. The light curve of D23H-1 shows evidence for additional variability in 2015, hence (Eq. 18) is not well-defined for this source. In the following section we therefore focus only on PTF-09ge and PTF-axc, since these two flares show a large flux increase with respect to a relatively well-defined baseline of late-time observations. We note that PTF-09ge appears to show a modest decline in flux at late times, suggesting our estimate of the difference flux of this flare could be slightly ( 10%) too low.
4. Parameter inference
We now discuss how to estimate the parameters of our dust reprocessing model (Sec. 2) from the observed light curve (Fig. 2). Our model assumes dust reprocessing in a thin shell, but this assumptions has little impact on our results (see Appendix). The results are summarized in Table 3.
To compare our model to the observations we use an equivalent form to Eq. 16,
The dust grains are not perfect black bodies and the IR spectrum of the dust is described by a modified Planck function:
with for grains of m and (Draine & Lee, 1984) and a constant. For convenience, we normalize this spectrum such that the amplitude of the reprocessing signal is dimensionless (i.e., ). The amplitude is determined from the observed IR luminosity; it measures the ratio between the observed IR luminosity and the amplitude of the reprocessing light curve, so we typically have (Fig. 3).
We use the well-sampled TDF PS1-10jh (Gezari et al., 2012) as a template for . Since we are in the regime where is much larger than the typical timescale of the optical flare (), the TDF energy determines the amplitude of the reprocessing signal, and the exact shape of the flare light curve is not important. The light curve of PS1-10jh provides a remarkably good match to PTF-09ge; only a small normalization shift (0.05 mag) and no temporal stretch have been applied. For the other flare, PTF-09axc, no late-time detections are available, and we will assume the post-flare light decays with the same rate as PTF-09ge.
Although our IR light curves (Fig. 2) are sparsely sampled, the amplitude and radius can be accurately determined as they are nearly orthogonal when is shorter than . In other words, our fiducial light curve model has very limited flexibility.
We use a maximum likelihood method to find the best-fit value of and . While we expect that the dust emits near the sublimation temperature ( K), lower temperatures are possible if dust only exists at radii larger than the sublimation radius. To account for the uncertainty due to the unknown dust temperature, we use both the W1 and W2 observations and we allow dust temperatures in the range . We maximize the logarithm of the likelihood given by:
with the observed IR luminosity and the sum running over the observations in the W1 and W2 band.
Contours of constant likelihood are shown in Fig. 3. We see that the shell radius of PTF-09ge can be accurately determined. For this flare, the probability at is vanishingly small, implying that the shape of the IR light curve is significantly different from the shape of the TDF light curve. For the other flare, PTF-09axc, an IR light curve that mirrors the TDF light curve can only be ruled-out at the 1 level.
To estimate the uncertainty on the relevant parameters, we sample the 3-dimensional likelihood using the straightforward Monte Carlo method of rejection sampling. We thus obtain distributions for sets of our three parameters (, , ) or for any scalar function that takes these parameters as input (i.e., and ). The 68% confidence interval for a given parameter is then given by the 16 and 84 percentiles of the distribution of the parameter.
Combinations of temperature and radius that yield (i.e., a total luminosity that is lower than the observed luminosity) are unphysical and excluded when computing the confidence intervals on and . The temperature range that is excluded by the requirement is below the lower limit on the temperature based on the W2 observations. Hence this cut mainly serves to the remove unphysical solutions near .
5.1. Origin of the IR emission
We can rule out that the observed transient IR emission originates from the same photosphere as the optical/UV emission of the tidal flare. First, the shape of the IR light curve is significantly different from the monotonic power-law decay that describes optical TDF emission. And second, the observed IR luminosity () is two orders of magnitude larger than the luminosity obtained from extrapolating the optical/UV black body spectrum, observed near the peak of the flare, to the time of the first WISE observation. In Fig. 4 we show the optical black body spectra of the two TDFs in our final sample. Since the optical observations are on the Rayleigh-Jeans side of the black body spectrum and well-fit with a single temperature, we can extrapolate this SED to find the expected WISE flux from the TDF photosphere at the time of the WISE observations. For both flares, this extrapolated flux is at least an order magnitude lower than the observed flux. We can thus conclude that the IR emission does not originate from the optical/UV photosphere. To match the observed W1 flux with a black body of K, the radius of this IR photosphere needs to be times larger than the radius of the optical/UV photosphere.
At the time of the first WISE observations, emission from hot dust exceeds the TDF emission for wavelengths nm (see Fig. 4) and therefore could be detectable with ground-based near-IR observations. For a generic dust reprocessing model, the IR emission is delayed, which could explain why for previous TDFs, a photosphere with K is not evident in optical spectra obtained near the peak of the flare.
Besides thermal emission from the TDF, a source of the observed transient IR flux could be synchrotron emission from a jet, such as observed for the relativistic TDF Swift J1644+57 (Burrows et al., 2011). However none of the TDFs in our sample have been detected in radio follow-up observations (van Velzen et al., 2013; Arcavi et al., 2014), ruling out emission from jets similar to Swift J1644+57.
With other potential sources of IR emission ruled out, dust reprocessing remains as the most plausible explanation for our observations. In fact, the light curve of PTF-09ge (Fig. 2) provides strong evidence for reprocessing by a spherical dust shell. As discussed in Section 2, when a spherical shell is briefly illuminated by a central point source, the area that is seen to emit simultaneously to an observer at large distance is constant with time.
After two observations with a near-constant flux, the third point in the 3.4 m light curve of PTF-09ge, 1.5 years after maximum light, shows a drop that is much steeper than the power-law decay of the optical TDF emission (Fig. 5). At this point, the peak of the flare emission, which lasted only a few months, has been reprocessed and all of the resulting IR photons have crossed the shell. The time delay between the TDF peak and this drop provides a direct measurement of the radius where the reprocessing happens. For a spherical shell, we find pc (Sec. 4).
The observed infrared luminosity of the two flares in our final sample is an order of magnitude fainter than recent theoretical estimates of dust emission from TDFs by Lu et al. (2016). This difference is partially due to the dust covering factor that Lu et al. adopted, which is an order of magnitude larger than the covering factor implied by our observations (see Sec. 5.3).
The light curve of the third flare with significant IR variability, D23H-1, cannot be explained by a single-shell reprocessing model because it shows a large change in flux at both 2 and 7 years after the optical peak. Future IR monitoring of this source is needed to measure the host galaxy baseline flux. We speculate that the IR light curve of D23H-1 can be explained by a concentration of dust at two different radii. Alternatively, the late-time re-brightening of D23H-1 could be explained if the source of reprocessed emission is the accretion disk of an AGN (which could be tested with additional optical monitoring of this source). Finally, we note that the late-time light curve of PTF-09ge shows evidence for low-level emission from regions beyond the sublimation radius of graphite grains.
5.2. TDF bolometric luminosity
With the radius of the dust shell measured, we can estimate the luminosity of the flare emission that has carved out this shell (, Eq. 12) in terms of the dust temperature. We anticipate the dust temperature is close to the sublimation temperature, K, but our observations at 4.6 m are not sensitive enough to provide a measurement of this temperature. This introduces an uncertainty on that is accounted for by marginalizing over temperature. As discussed in detail in Section 4, we use the observations at 3.4 and 4.6 m to find the likelihood of the model light curve as a function of and . For the flare PTF-09ge, we find (68% confidence interval).
The TDF luminosity inferred from a dust reprocessing light curve closely approximates the bolometric flare luminosity because dust absorption is efficient from optical to soft X-ray frequencies, thus covering the full spectral range where thermal TDFs emit most of their energy — including the extreme-UV which is nearly impossible to observe directly due to absorption by neutral hydrogen. The ratio between the luminosity inferred from the reprocessing light curve and the black body luminosity estimated from the optical observations is for PTF-09ge (logarithmic units). Either of two plausible explanations may account for this large bolometric correction222For AGN, the bolometric correction refers to the ratio between at a given frequency and the total luminosity. Since optical/UV observations of TDF are well-described by a black body spectrum, it is more instructive to define the TDF bolometric correction with respect to the black body luminosity..
First, the bolometric correction may be due to dust extinction in the host galaxy (i.e., due to dust along our line of sight to the black hole). This reddening decreases the observed black body temperature and, since the flare SED peaks at UV wavelengths, significantly reduces the inferred black body luminosity. The observed slope of the optical SED of PTF-09ge () is more shallow than the slope for the Rayleigh-Jeans limit, hence the intrinsic temperature of this flare could be higher than the observed temperature. Using the Calzetti et al. (2000) extinction law to correct the observed optical spectrum of PTF-09ge for a dust column with increases the black body luminosity by a factor . An analysis of sodium absorption in the spectra of type 1 AGN suggests that an extinction of to a galactic nucleus is not uncommon (Baron et al., 2016).
Alternatively, if the extinction to the center of the galaxy is small, the bolometric correction could be explained by adding an X-ray emitting component to the flare SED. The TDF ASASSN-14li (Holoien et al., 2015; Miller et al., 2015; van Velzen et al., 2016) showed two black body spectra, with temperatures K and keV. The latter dominates the total energy output and could thus explain the bolometric correction to the optical luminosity of PTF-09ge.
If the accretion of stellar debris is radiatively efficient (), our estimate of the total radiated energy of PTF-09ge implies an upper bound on the accreted mass of . Recent numerical simulations (Shiokawa et al., 2015) show a similar mass accretion after the disruption of a solar type star, hence our observations are consistent with high radiative efficiency in super-Eddington accretion disks (Jiang et al., 2014). This inference also points to a full disruption of a solar-type star, while a partial disruption would be inferred if no bolometric correction is applied to the light curve.
Note. – For the two tidal flares with transient 3.4 m flux, we list the observed peak IR luminosity in the second column. The third to last columns show the results of applying our dust reprocessing model to the IR light curve. From the break in the light curve we find the radius of the reprocessing shell (). For a given dust temperature, we can estimate the peak luminosity of the tidal flare over the frequencies where dust absorbs the flare photons (). Integrating over the flare light curve yields the total radiated energy (). Finally, the covering factor () is given by the total energy radiated by the dust divided by . We list the parameters that correspond to the maximum likelihood and the uncertainties reflect the 68% confidence interval. This interval includes a marginalization over the dust temperature, which provides the largest contributions to the uncertainty.
5.3. Dust covering factor
Our observations are the first to probe dust within 0.1 parsec of the center of non-active galaxies. We can use the total energy radiated in the IR () and the energy able to heat graphite dust () to find the covering factor of this dust, . For both TDFs in our final sample we find % (Table 3).
The fractional uncertainty on is smaller than for since the former has a weaker temperature-dependence. The ratio between the total IR luminosity and the portion we observe at 3.4m (cf. Eq. 5), increases with temperature, for temperatures K, but somewhat less steeply at lower temperatures. Since (Eq. 12), the fractional uncertainty on the covering factor is a factor smaller than the fractional uncertainty on .
Galaxy-to-galaxy fluctuations in the dust size distribution are expected to have a small influence on the inferred covering factor. Dust lanes in E/S0 galaxies have an extinction curve that is similar to the Milky Way (Finkelman et al., 2012), which implies a similar peak of the grain size distribution (Goudfrooij et al., 1994). For a sample of 26 early-type galaxies, the mean grain size difference with respect to the Milky Way is 8% (Patil et al., 2007). Since (Eq. 12), this fluctuation of the grain size translates to an uncertainty of 0.1 dex on .
While the dust distribution at the centers of galaxies is not constrained by observations, most mechanisms that can alter the distribution (e.g., sputtering) will act to reduce the number of small grains relative to large grains and therefore not affect our estimate of the covering factor, unless there are also agglomerative mechanisms particular to galactic nuclei.
The information contained in the IR light curve (Fig. 5) is not sufficient to constrain to what degree the dust geometry departs from spherical symmetry (see Appendix). However this does not affect our ability to measure the covering factor because this parameter is a measure of the absorbed energy, which is independent of the dust geometry (this is demonstrated in Fig. 6, right panel).
The covering factor inferred from our observations is almost two orders of magnitude smaller than typical dusty tori in Seyfert galaxies (Barvainis, 1987). This is not surprising since a high accretion rate is likely required to build and sustain a torus that covers a large solid angle (e.g., Pier & Krolik, 1992), while the host galaxies of the TDFs show no signs of high accretion rates prior to the stellar disruption. The Galactic Center likely provides a better comparison. Interestingly, the nuclear dust in the TDF host galaxies is different from the circumnuclear ring of molecular gas at the Galactic Center; this ring has a covering factor of about 20% and a sharp inner edge at pc from Sgr A* (Genzel et al., 2010; Ferrière, 2012). Inside this edge, a region known as the “ionized cavity”, free floating dust particles will be sublimated by the UV radiation of the nuclear star cluster. Yet dust can exist within dense molecular cores that are observed inside the cavity (Christopher et al., 2005). Similar clumps may be the source of nuclear dust in the TDF host galaxies. Alternatively, if the TDF host galaxies have a nuclear star cluster with an old stellar population or lack a nuclear star cluster, dust on 0.1 pc scales may simply originate from streams of cold gas that fall toward the central black hole without being ionized.
Finally, we point to a potential section effect which could explain the relatively small covering factor of the galaxies in our sample. Since the optical/UV SED of the TDF candidates found so far is relatively blue ( K), extinction along the line of sight will quickly reduce their detectability. While color is not used as a selection criterion to find TDF candidates in optical transient surveys (e.g., van Velzen et al., 2011; Arcavi et al., 2014), these surveys might still have a moderate bias to finding flares in galaxies with a small amount of dust. If the amount of dust along our line of sight to the galaxy center is strongly correlated with the amount of circumnuclear dust on sub-parsec scales, a TDF sample could yield a biased view of the nuclear dust covering factor. For a large sample of TDFs, this potential bias could be quantified by comparing the dust covering factor, as measured from the reprocessing light curve, to other dust extinction estimates (e.g., as obtained from the Balmer decrement or stellar population synthesis).
6. Conclusions and Outlook
Our main conclusions are:
We have discovered variable 3.4 m emission for three of the five tidal disruption flares in our sample (Fig. 5), the host-subtracted IR luminosity is .
For two of the five TDFs in our sample, the observed IR light curves can be explained by reprocessing of UV/X-ray emission from the flare in a thin shell at pc from the black hole (Fig. 5).
The radius of the IR-emitting region yields an estimate of the bolometric luminosity of the TDF. We find a typical luminosity of erg s (Table 3), which is a factor larger than the observed black body luminosity.
The covering factor of the nuclear dust that is responsible for reprocessing the TDF light is % (Table 3).
For a dust temperature of K, ground-based IR observations of newly discovered TDFs could be used to estimate the dust temperature and obtain reprocessing light curves at higher cadence than presented in this work. These observations would reduce the uncertainty on the TDF bolometric luminosity and will also provide information on the geometry of the nuclear dust. High-resolution spectroscopic observations (e.g., with JWST) can be used to study in detail the response of dust particles to a sudden burst of UV radiation.
- Arcavi et al. (2014) Arcavi, I., et al. 2014, ApJ, 793, 38
- Asplund et al. (2009) Asplund, M., Grevesse, N., Sauval, A. J., & Scott, P. 2009, ARA&A, 47, 481
- Baron et al. (2016) Baron, D., Stern, J., Poznanski, D., & Netzer, H. 2016, ArXiv e-prints: 1603.06948
- Barvainis (1987) Barvainis, R. 1987, ApJ, 320, 537
- Berger et al. (2012) Berger, E., Zauderer, A., Pooley, G. G., Soderberg, A. M., Sari, R., Brunthaler, A., & Bietenholz, M. F. 2012, ApJ, 748, 36
- Bloom et al. (2011) Bloom, J. S., et al. 2011, Science, 333, 203
- Burrows et al. (2011) Burrows, D. N., et al. 2011, Nature, 476, 421
- Calzetti et al. (2000) Calzetti, D., Armus, L., Bohlin, R. C., Kinney, A. L., Koornneef, J., & Storchi-Bergmann, T. 2000, ApJ, 533, 682
- Cenko et al. (2012) Cenko, S. B., et al. 2012, ApJ, 753, 77
- Chornock et al. (2014) Chornock, R., et al. 2014, ApJ, 780, 44
- Christopher et al. (2005) Christopher, M. H., Scoville, N. Z., Stolovy, S. R., & Yun, M. S. 2005, ApJ, 622, 346
- Dou et al. (2016) Dou, L., Wang, T.-g., Jiang, N., Yang, C., Lyu, J., & Zhou, H. 2016, ArXiv e-prints
- Draine & Lee (1984) Draine, B. T., & Lee, H. M. 1984, ApJ, 285, 89
- Esquej et al. (2008) Esquej, P., et al. 2008, A&A, 489, 543
- Ferrière (2012) Ferrière, K. 2012, A&A, 540, A50
- Finkelman et al. (2012) Finkelman, I., Brosch, N., Funes, J. G., Barway, S., Kniazev, A., & Väisänen, P. 2012, MNRAS, 422, 1384
- Genzel et al. (2010) Genzel, R., Eisenhauer, F., & Gillessen, S. 2010, Reviews of Modern Physics, 82, 3121
- Gezari et al. (2006) Gezari, S., et al. 2006, ApJ, 653, L25
- Gezari et al. (2009) —. 2009, ApJ, 698, 1367
- Gezari et al. (2012) —. 2012, Nature, 485, 217
- Gezari et al. (2013) —. 2013, ApJ, 766, 60
- Goudfrooij et al. (1994) Goudfrooij, P., de Jong, T., Hansen, L., & Norgaard-Nielsen, H. U. 1994, MNRAS, 271, 833
- Guhathakurta & Draine (1989) Guhathakurta, P., & Draine, B. T. 1989, ApJ, 345, 230
- Holoien et al. (2014) Holoien, T. W.-S., et al. 2014, MNRAS, 445, 3263
- Holoien et al. (2015) —. 2015, MNRAS, in press
- Jiang et al. (2016) Jiang, N., Dou, L., Wang, T., Yang, C., Lyu, J., & Zhou, H. 2016, ArXiv e-prints
- Jiang et al. (2014) Jiang, Y.-F., Stone, J. M., & Davis, S. W. 2014, ApJ, 796, 106
- Jones et al. (2015) Jones, D. O., Scolnic, D. M., & Rodney, S. A. 2015, PythonPhot: Simple DAOPHOT-type photometry in Python, Astrophysics Source Code Library
- Kaspi et al. (2000) Kaspi, S., Smith, P. S., Netzer, H., Maoz, D., Jannuzi, B. T., & Giveon, U. 2000, ApJ, 533, 631
- Komossa & Bade (1999) Komossa, S., & Bade, N. 1999, A&A, 343, 775
- Komossa et al. (2008) Komossa, S., et al. 2008, ApJ, 678, L13
- Koshida et al. (2014) Koshida, S., et al. 2014, ApJ, 788, 159
- Lang (2014) Lang, D. 2014, AJ, 147, 108
- Lang et al. (2014) Lang, D., Hogg, D. W., & Schlegel, D. J. 2014, ArXiv e-prints: 1410.7397
- Law et al. (2009) Law, N. M., et al. 2009, PASP, 121, 1395
- Levan et al. (2011) Levan, A. J., et al. 2011, Science, 333, 199
- Levan et al. (2016) —. 2016, ApJ, 819, 51
- Lu et al. (2016) Lu, W., Kumar, P., & Evans, N. J. 2016, MNRAS, 458, 575
- Mainzer et al. (2011) Mainzer, A., et al. 2011, ApJ, 731, 53
- Mainzer et al. (2014) —. 2014, ApJ, 792, 30
- Martin et al. (2005) Martin, D. C., et al. 2005, ApJ, 619, L1
- Masci (2013) Masci, F. 2013, ArXiv e-prints: 1301.2718
- Masci & Fowler (2009) Masci, F. J., & Fowler, J. W. 2009, in Astronomical Society of the Pacific Conference Series, Vol. 411, Astronomical Data Analysis Software and Systems XVIII, ed. D. A. Bohlender, D. Durand, & P. Dowler, 67
- Miller et al. (2015) Miller, J. M., et al. 2015, Nature, 526, 542
- Palaversa et al. (2016) Palaversa, L., Gezari, S., Sesar, B., Stuart, J. S., Wozniak, P., Holl, B., & Ivezić, Ž. 2016, ApJ, 819, 151
- Patil et al. (2007) Patil, M. K., Pandey, S. K., Sahu, D. K., & Kembhavi, A. 2007, A&A, 461, 103
- Peng et al. (2002) Peng, C. Y., Ho, L. C., Impey, C. D., & Rix, H.-W. 2002, AJ, 124, 266
- Pier & Krolik (1992) Pier, E. A., & Krolik, J. H. 1992, ApJ, 399, L23
- Roth et al. (2015) Roth, N., Kasen, D., Guillochon, J., & Ramirez-Ruiz, E. 2015, ArXiv e-prints: 1510.08454
- Saxton et al. (2012) Saxton, R. D., Read, A. M., Esquej, P., Komossa, S., Dougherty, S., Rodriguez-Pascual, P., & Barrado, D. 2012, A&A, 541, A106
- Shiokawa et al. (2015) Shiokawa, H., Krolik, J. H., Cheng, R. M., Piran, T., & Noble, S. C. 2015, ApJ, 804, 85
- Stoughton et al. (2002) Stoughton, C., et al. 2002, AJ, 123, 485
- van Velzen et al. (2013) van Velzen, S., Frail, D. A., Körding, E., & Falcke, H. 2013, A&A, 552, A5
- van Velzen et al. (2015) van Velzen, S., Gorjian, V., Krolik, J., & Mendez, A. 2015, Dust Reprocessing of Stellar Tidal Disruption Flares, Spitzer Proposal
- van Velzen et al. (2011) van Velzen, S., et al. 2011, ApJ, 741, 73
- van Velzen et al. (2016) —. 2016, Science, 351, 62
- Waxman & Draine (2000) Waxman, E., & Draine, B. T. 2000, ApJ, 537, 796
- Weingartner & Draine (2001) Weingartner, J. C., & Draine, B. T. 2001, ApJ, 548, 296
- Wright et al. (2010) Wright, E. L., et al. 2010, AJ, 140, 1868
- Yoon et al. (2015) Yoon, Y., et al. 2015, ApJ, 808, 96
- York et al. (2000) York, D. G., et al. 2000, AJ, 120, 1579
- Zauderer et al. (2011) Zauderer, B. A., et al. 2011, Nature, 476, 425
Appendix A Alternative dust geometries
Our fiducial model (Eq. 4) assumes emission from a single shell. In this section we explore to what extent deviations from this assumption are constrained by the data and how they affect the inferred parameters.
We first consider a superposition of shells at different radii. The reprocessing light curve for this geometry follows by integrating Eq. 4 over ,
Here is the number density of the dust as a function of radius and , with the temperature of the inner shell (cf. Eq. 12). Eq. A1 does not account for self-shielding of the dust, which is justified as long as the covering factor is small. The width of the response function of each shell (Eq. 16) is proportional to the shell radius, hence the amplitude of the reprocessing light curve decreases with , multiplied by a small spectral correction of order . As a result, the IR light curve for a constant dust density has a much more shallow flux decay at compared to a single shell. The light curve of PTF-09ge appears inconsistent with a constant density (Fig. 6, left panel). Using only the W1 observations, the log likelihood difference compared to the fiducial model is .
Clearly, our observations have too few degrees of freedom to measure a slope of a dust density profile, , but we can conclude that steep profiles () provide a better description of the data for PTF-09ge. For example, a Bondi density profile () with an inner radius at pc yields a reasonable fit (log likelihood difference with respect to the fiducial model is -1.4).
Finally, we also consider a disk geometry for the dust distribution. For a thin disk, the reprocessing light curve is given by
Here parametrizes the coordinates of disk, which run between , with the inclination of the disk (e.g., for a face-on disk, , all reprocessed emission has the same time delay of ). Our observations cannot disentangle the degeneracy between radius and inclination of a disk geometry. Yet we can use Eq. A2 to fit for the radius of the disk for a range of inclinations and then use Eq. 12 to estimate the absorbed luminosity and the dust covering factor for these radii. We find that the range of inclinations that are consistent with the data of PTF-09ge yield a similar value of the covering factor compared to our fiducial dust model (Fig. 6, right panel). This should not be a surprise since the covering factor measures an energy ratio, i.e., it has been integrated over the geometry of the emitting region.