X-RAY AND OPTICAL CORRELATION OF TYPE I SEYFERT NGC 3516 STUDIED WITH SUZAKU AND JAPANESE GROUND-BASED TELESCOPES

The X-ray monitoring observations of NGC 3516 were performed in the Suzaku AO-8 cycle during 2013–2014. Specifically, NGC 3516 was observed with Suzaku five times from 2013 April 10 to May 29, with intervals of ∼1–2 weeks (see Table 1). Another pointing was carried out on November 7, and the last one on 2014 April 4, making a total of seven observations, with various intervals from days to months. In the present paper, we call the first, second, ... and seventh observations epoch 1, epoch 2, ... and epoch 7, respectively. The exposure in each epoch was ∼50 ks, except for epoch 2 which had an exposure of ∼20 ks. In all these epochs, the XIS (Koyama et al. 2007) and HXD (Takahashi et al. 2007) were operated in their normal modes, and data from the XIS and HXD-PIN are utilized in the present paper, except for epoch 7, in which the HXD-PIN count rate was too low to significantly detect the signals.

The data of the XIS and HXD-PIN were processed via the software version 2.4. In the XIS analysis, we added the data of XIS0 and XIS3, and utilized these as FI data, while we did not use XIS1 due to their higher background. The XIS source events were extracted from a circle of radius 120'' with its center on the source. Background events were accumulated from a surrounding annular region with an inner radius of 180'' and outer radius of 270''. The response and ancillary-response files were prepared via software in HEASOFT 6.14 called xisrmfgen and xissimarfgen (Ishisaki et al. 2007), respectively. The HXD data were prepared in the same way as those from the XIS. In the HXD analysis, non-X-ray background and cosmic X-ray background were estimated using standard models described by Fukazawa et al. (2009) and Gruber et al. (1999), respectively. They were subtracted from the on-source data.

The optical photometric monitoring observations of NGC 3516 were performed by using the charge-coupled device (CCD) cameras installed at the Pirka, MITSuME, Kiso Schmidt, Nayuta, and Kanata telescopes, of which the major parameters are summarized in Table 2. Simultaneous monitoring observations were performed at the same seven epochs as the X-ray observations by Suzaku. In addition, from 2013 January to 2014 April, these telescopes densely monitored the source even without simultaneous X-ray coverage, with a typical intervals of ∼1 day. In the present paper, we present the B-band and g'-band photometric data, because the optical continuum emission of AGNs is generally bluer than that of the host galaxy, thus ensuring better sensitivity to the AGN signals. The images were reduced using IRAF¹⁶ following the standard procedures of image reduction for CCD detectors.

Figure 1 shows light curves in the 2–3 and 3–10 keV bands in all seven epochs. While we can see gradual flux changes, the source did not show short-term variations, on timescales of several hours, that were observed by XMM-Newton in 2006 (Mehdipour et al. 2010). The 2–3 keV light curves are clearly synchronized with those in the 3–10 keV band. The highest flux was obtained in epoch 4 and the lowest in epoch 7, with the amplitude of peak-to-peak variation reaching a factor of ∼20. To quantify the flux variability, we made, in Figure 2, a count–count plot between count rates in the 2–3 and 3–3.5 keV bands, where these objects generally exhibit rather high variability (Noda et al. 2013a). The count–count plot of NGC 3516 derived in 2009 October (Noda et al. 2013b) is also shown for reference, after correcting the data for a difference in pointing position and long-term changes in the detector response. Surprisingly, all the data points in 2013–2014 line up with those in 2009, thus defining an almost linear correlation between the two bands. This means that the variable component in those energy bands observed in 2013–2014 has nearly the same spectral shape as that detected in 2009.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In energy bands higher than 3.5 keV, a reflection component with a prominent Fe–Kα emission line at ∼6.4 keV is expected to become more dominant than in energies below 3.5 keV. Therefore, this variable component should be separated from the reflection via spectral fitting. Although the reflection component is considered to be almost stationary within a week (Noda et al. 2013a), it can vary on timescales of months, because a broad-line region and/or a dusty torus, where the reflection emission possibly originates, are known to be located at several tens to hundreds of light-days (e.g., Koshida et al. 2014). Accordingly, we treat the reflection signal as a variable component.

As presented in Figure 3, we produced seven time-averaged and background-subtracted spectra, one from each epoch, and performed a simultaneous model fitting to them. Because of the presence of diffuse X-ray emissions in the host galaxy (e.g., Costantini et al. 2000; George et al. 2002) and the effects of variations in absorption (e.g., Turner et al. 2011), we ignored energy bands lower than 2 keV, and utilized the 2–45 keV energies in each spectrum. The employed spectral model is wabs ∗ (cutoffpl + pexmon) in XSPEC12, which is consistent with that utilized in Noda et al. (2013a, 2013b). Here, wabs is a model for the photoelectric absorption (Morrison & McCammon 1983), and its column density N_H was allowed to differ among the epochs. The cutoffpl model represents the PL continuum with a high-energy exponential rolloff, utilized to emulate the inverse Comptonization radiation. The reflection continuum and the Fe–Kα emission line, both produced by this continuum, are reproduced by pexmon (Nandra et al. 2007).

Zoom In Zoom Out Reset image size

In the fitting, the photon index Γ was tied among the epochs (because of Figure 2) and left free. The cutoff energy E_cut was fixed at 150 keV based on the typical value reported by Malizia et al. (2014), while the normalization N_PL was left free in each epoch. For pexmon, Γ and E_cut of the incident PL were tied to those in cutoffpl, and the Fe abundance A_Fe, inclination angle I, and reflection fraction f_ref were fixed at 1 Solar, 60°, and 1, respectively. The normalization N_ref of the pexmon component was allowed to differ among the epochs. The N_ref value is determined almost solely by the Fe–Kα line intensity, because the ratio between the Fe–Kα line intensity and the reflection continuum in pexmon does not change when E_cut, A_Fe, and I are all fixed. As summarized in Figure 3 and Table 3, the simultaneous fitting was successful with χ²/dof = 1184.0/1113. As expected from Figure 2, all the spectra have been reproduced by a hard PL continuum having Γ ∼ 1.75, which is consistent with that in 2009 (${\rm{\Gamma }}={1.72}_{-0.12}^{+0.08}$), together with a reflection component accompanied by a prominent Fe–Kα emission line. Although N_H varied slightly, the low-energy shapes of the time-averaged spectra were not so strongly affected (see Figure 3).

Table 3.  Parameters Obtained by Simultaneous Fitting to all the 2–45 keV Time-averaged Spectra

Component	Parameter	Epoch 1	Epoch 2	Epoch 3	Epoch 4	Epoch 5	Epoch 6	Epoch 7
wabs	NHa	${1.93}_{-0.15}^{+0.14}$	${1.86}_{-0.11}^{+0.10}$	${1.21}_{-0.07}^{+0.06}$	${0.94}_{-0.06}^{+0.05}$	1.16 ± 0.07	${1.52}_{-0.17}^{+0.16}$	${2.10}_{-0.25}^{+0.21}$
cutoffpl	ΓPL	${1.75}_{-0.02}^{+0.01}$
	${E}_{\text{cut,PL}}$ (keV)	150 (fix)
	NPLb	${1.00}_{-0.03}^{+0.02}$	${3.94}_{-0.13}^{+0.11}$	${5.44}_{-0.12}^{+0.09}$	${6.93}_{-0.15}^{+0.12}$	${4.10}_{-0.09}^{+0.07}$	${0.69}_{-0.03}^{+0.02}$	0.39 ± 0.02
	FPLc	0.37 ± 0.01	${1.44}_{-0.05}^{+0.04}$	${1.99}_{-0.04}^{+0.03}$	${2.53}_{-0.05}^{+0.04}$	${1.50}_{-0.03}^{+0.02}$	0.25 ± 0.01	0.14 ± 0.01
pexmon	Γref	= ΓPL
	${E}_{\text{cut,ref}}$ (keV)	= ${E}_{\text{cut, PL}}$
	fref	1 (fix)
	AFe (Z⊙)	1 (fix)
	Nrefd	${3.88}_{-0.26}^{+0.25}$	${4.24}_{-0.50}^{+0.51}$	${5.17}_{-0.44}^{+0.42}$	${5.88}_{-0.49}^{+0.47}$	${5.60}_{-0.41}^{+0.40}$	2.73 ± 0.21	2.51 ± 0.18
	Frefe	1.10 ± 0.07	${1.20}_{-0.14}^{+0.15}$	1.46 ± 0.12	${1.67}_{-0.14}^{+0.13}$	${1.59}_{-0.12}^{+0.11}$	0.78 ± 0.06	0.71 ± 0.05
	Ftotalf	0.48 ± 0.04	1.56 ± 0.07	${2.14}_{-0.07}^{+0.06}$	2.70 ± 0.07	${1.66}_{-0.06}^{+0.05}$	0.33 ± 0.03	0.21 ± 0.02
	χ2/dof	1184.0/1113

Notes. Errors refer to 1σ confidence ranges.

^aEquivalent hydrogen column density in units of 10²² cm⁻². ^bThe power-law normalization at 1 keV, in units of 10⁻³ photons keV⁻¹ s⁻¹ cm⁻² at 1 keV. ^cThe 2–10 keV flux of the PL component without being absorbed, in units of 10⁻¹¹ erg s⁻¹ cm⁻². ^dThe pexmon normalization at 1 keV, in units of 10⁻³ photons keV⁻¹ s⁻¹ cm⁻² at 1 keV. ^eThe 2–10 keV flux of the reflection component without being absorbed, in units of 10⁻¹² erg s⁻¹ cm⁻². ^fThe 2–10 keV total flux without being absorbed in units of 10⁻¹¹ erg s⁻¹ cm⁻², calculated as a sum of F_PL and F_ref.

Download table as:  ASCIITypeset image

The highest and lowest 2–10 keV fluxes were ∼2.7 ×10⁻¹¹ erg s⁻¹ cm⁻² in epoch 4 and ∼0.3 × 10⁻¹¹ erg s⁻¹ cm⁻² in epoch 7, respectively. They are in a range of just 7%–70% of the 2–10 keV flux of ∼4.0 × 10⁻¹¹ erg s⁻¹ cm⁻² averaged over 1997–2002 (RXTE AGN Timing & Spectral Database; Maoz et al. 2002); thus NGC 3516 was in an X-ray-faint phase during the present observations. Table 3 further gives the 2–10 keV fluxes of the hard PL continuum and the reflection components, calculated separately, after removing the absorption. The 2–10 keV light curves of the two spectral components, derived in this way, are presented in Figure 4(top). The time of each data point refers to the middle epoch (in MJD) of that observation as given in Table 1. Interestingly, the flux of the reflection component changed significantly by a factor of 2 on a timescale of several months. However, the hard PL component varied by almost an order of magnitude or more in amplitude. Clearly, the large change in intensity in Figure 2 was mostly carried by the hard PL variation. The largest flux change in the hard PL emission in 2–10 keV is a difference of ∼2.4 × 10⁻¹¹ erg s⁻¹ cm⁻², and a factor of ∼18 between epochs 4 and 7.

Zoom In Zoom Out Reset image size

To get further information about the AGN activity in NGC 3516, we also performed the spectral analysis including the energy band below 2 keV. For this purpose we selected epoch 7, where the AGN was faintest. In fact, as shown in Figure 5, the spectrum on this occasion exhibits a prominent soft X-ray excess at ≲2 keV, which is much less prominent in the other epochs. In order to examine whether the soft excess structure originates from the AGN activity, we fitted the 0.5–10 keV spectrum at epoch 7 by two contrasting models: (a) wabs ∗ (cutoffpl + pexmon) + apec and (b) wabs ∗ (cutoffpl + pexmon + kdblur ∗ reflionx). The soft excess is explained in model (a) by apec (Smith et al. 2001), which represents thin thermal plasma emission from the host galaxy, and in model (b) by a kdblur ∗ reflionx model (Laor 1991; Ross & Fabian 2005), which describes a relativistically smeared reflection continuum. As a result, the fit with model (a) was successful with χ²/dof = 100.3/84 (Figure 5(a)), in which the soft excess, apparently involving emission lines, is reproduced by the plasma emission model. On the other hand, the fit with model (b) was unsuccessful with χ²/dof = 355.9/82, mainly due to the lack of the emission lines at ∼0.85 keV and a convex data shape in the 2–5 keV band (Figure 5(b)). Therefore, the spectrum below 2 keV is dominated by the host galaxy emission, at least in epoch 7. This justifies us in limiting the spectral studies to energies above ∼2 keV. Of course, the AGN emission, when bright, probably dominates down to ∼0.5 keV, but it is strongly affected by ionized-absorption features at ∼0.8 keV, as shown in Figure 5(a).

Zoom In Zoom Out Reset image size

According to changes in atmospheric seeing during an observation, the flux of an AGN within an aperture changes differently from that of the host galaxy, leading to a significant uncertainty in measuring the flux variation of a faint AGN hosted by a bright galaxy. To minimize such uncertainty in the photometry of the nucleus of NGC 3516, we therefore performed Difference Image Photometry (DIP: Crotts 1992; Tomaney & Crotts 1996). As presented in Figure 6, for example, a reference image was created by stacking many images obtained on the same night with good seeing conditions. Then, we matched the position, the photometric intensity, and the point-spread function (PSF) of the reference image to those of each individual image and subtracted the former from the latter. Two nearby field stars, located at (R.A., decl.) = (11:06:28.51, +72:35:46.2) and (11:07:00.41, +72:33:33.3), were used as the reference to match the photometric intensity and the PSF. The IRAF psfmatch task was used for matching the PSFs. After that, for each individual image, the residual flux at the center of the galaxy was measured relative to the two nearby reference stars, with a circular aperture of ϕ = 4 × FWHM of the PSF in diameter. Finally, the residual flux data on the same observing night were averaged to obtain the flux difference of the nucleus of NGC 3516 at that epoch with respect to the reference image. The photometric error was estimated from the ensemble scatter of the DIP fluxes and the number of the images obtained on the same night. These procedures were applied to the data set from each telescope individually.

Zoom In Zoom Out Reset image size

Because the observing epochs of the reference images used in the DIP analysis depend on the telescopes, there can be systematic offsets between them. Systematic differences in the scaling factor among them are also possible because of the differences in the filter color term. To match the offset and the scaling factor among the four telescopes, we therefore performed a linear regression analysis on their B-band flux data sets from the same observing nights. An example, between two telescopes, is presented in Figure 7. The reduced χ² of the linear regression was much larger than unity when only the photometric errors were incorporated. This suggests the presence of some systematic errors, so we added an error σ_add to the photometric errors in quadrature, and regarded it as the systematic error of the photometry. By requiring the reduced χ² to become unity, we obtained σ_add ∼ 0.06 mJy.

Zoom In Zoom Out Reset image size

The B- and V-band magnitudes of the two reference stars were calibrated relative to those of the more distant field stars whose magnitudes are presented in Sakata et al. (2010); these are based on the wide-field image data of the KWFC, in which both stars were observed simultaneously. The g'-band magnitudes of the reference stars were estimated from their B- and V-band magnitudes according to Jordi et al. (2006).

The light curves of the nucleus of NGC 3516, thus derived in the optical B and g' bands, are presented in the bottom panel of Figure 4, and the data are listed in Table 4, without correction for the Galactic extinction. They are both the differential fluxes with respect to the reference image obtained from the image data on 2014 April 8 (UTC), and the error bars of the data points do not include σ_add. The total numbers of photometric data points are 180 and six for the B and g' bands, respectively. As shown in Figure 4, the B-band flux varied rather similarly to that of the hard PL component in the 2–10 keV band. Although the number of data points is small, the variation in g'-band flux also followed them. The larger relative scatters of the g'-band data points are caused by the small telescope size and the large pixel scale of the camera (MITSuME) with which the data were obtained. These results indicate that the optical continuum and the 2–10 keV PL component varied in a good correlation with each other. The B-band flux density varied by ∼1.2 mJy between the peak and bottom. After correcting for the Galactic extinction and the optical extinction of the nucleus of NGC 3516 (A_B ∼ 1.3 mag in total), this becomes ∼4.0 mJy, corresponding to ∼2.7 × 10⁻¹¹ erg s⁻¹ cm⁻² in νF_ν units.¹⁷

To complete this DIP analysis, we estimated the AGN flux at the faintest phase during our observation. We applied an aperture photometry to the stacked B-band image obtained on 2014 April 7 by the MINT attached on the Nayuta telescope (Figure 6), because on this day it achieved the best seeing among the four. Sakata et al. (2010) estimated the host-galaxy flux in the B band as 8.38 ± 0.18 mJy,¹⁸ within an aperture diameter of ϕ = 8.3 arcsec with the sky reference of a ϕ = 11.1–13.9 arcsec annulus. By applying the photometry with the same parameters, we obtained the B-band flux density of 8.44 mJy after the correction for the Galactic extinction, yielding an AGN flux of 0.1 ± 0.2 mJy. Considering the various systematic errors in the photometry, the B-band AGN flux at the faintest phase during the observations is estimated to be about a few times 0.1 mJy, with a similar amount of flux error.

As shown in Figure 4, the optical flux changed almost simultaneously with the 2–10 keV hard PL flux. In order to examine possible time lags between them, we first focus on their expanded light curves from epochs 1 to 5, as presented in Figure 8(a). Clearly, the AGN became brightest in both bands at epoch 4. Also plotted in the same figure is the X-ray light curve shifted by +2 days and +4 days. Apparently, the X-ray light curve with a shift of a few days shows the closest agreement with the B-band light curve, suggesting that the variation of the X-rays preceded that of the B band by a few days.

Zoom In Zoom Out Reset image size

Figure 8(b) presents the correlation of the X-ray flux with temporal shifts of 0, +2, and +4 days, against the B-band flux at the delayed epochs. The latter was estimated by averaging the B-band fluxes in close observations within ±1 day. The correlation between those fluxes appears to be strongest when the X-ray light curve is shifted by +2 days. This reconfirms that the X-ray variations preceded those in the B band by a few days. Below, let us quantify the suggested time lag with two methods.

3.3.1. Cross-correlation Analysis

One method is the interpolated cross-correlation function (ICCF) method, which has been widely used (White & Peterson 1994; Peterson et al. 1998). This is the same as the ordinary cross-correlation function (CCF), but either or both sets of time-series data are interpolated, so that the CCF can be calculated even if the two data sets have different samplings, or when either or both are irregularly sampled.

Since the monitoring cadence of the B-band data was much higher than that of the X-ray data, only the B-band light curve was interpolated, to make it nearly continuous. The interpolated B-band light curve was calculated using a fitting code developed by Zu et al. (2011), which assumes a damped random walk (DRW) model for the flux variation. The DRW model has been demonstrated to be a good statistical model of flux variations of the UV–optical continuum emission of AGNs (e.g., Kelly et al. 2009; Kozłowski et al. 2010; MacLeod et al. 2010, 2012; Zu et al. 2013). In this interpolation, the systematic error σ_add of the photometry, either as determined in Section 3.2 or somewhat changed, was added in quadrature.

Figure 9 presents the observed B-band flux data and their interpolation. Thus, the interpolated light curve does not follow the observed data when we employ σ_add = 0.06 mJy as determined in Section 3.2: the peak flux of the interpolated light curve was lower and its decrease in flux after the peak was slower. We generally found that a larger value of σ_add made the interpolated light curve smoother and less variable, probably because a larger σ_add would work as if applying a stronger low-pass filter. Therefore, we performed the CCF analysis using the B-band data with σ_add = 0 mJy in addition to 0.06 mJy to examine the possible uncertainty in lag caused by σ_add. As shown in Figure 9, the interpolated light curves then follow the observed data much better when σ_add is reduced.

Zoom In Zoom Out Reset image size

The X-ray data of epochs 1–5 and the optical data of MJD = 56316–56511 were selected for the CCF analysis, because the most remarkable flux variations were present in those epochs and also because the B-band light curve data were sampled densely (Figures 4 and 8). The calculated CCFs, in which the B-band data were smoothed using σ_add = 0 and 0.06 mJy, are presented in Figure 10, in which the time lag τ is defined to be positive if the optical flux precedes the X-rays. As shown in Figure 10, the CCF is peaked at about τ = −2 days. This reconfirms the inference from Figure 8 and implies that the optical variation is delayed from that in X-rays by ∼2 days. The CCF value at the peak is very high, >0.98, as indicated by Figure 8(b).

Zoom In Zoom Out Reset image size

For quantitative measurement of τ, we use the centroid of the CCF peak, τ_cent, which is computed from all neighboring points around the CCF peak where CCF values are >0.95 times that of the peak. The uncertainty of τ_cent is estimated using the model-independent Monte Carlo method of flux randomization and random subset sampling (FR/RSS) introduced by Peterson et al. (1998, 2004). The FR method modifies the observed fluxes in each realization randomly within the errors assigned to the individual data points, and the RSS method randomly extracts the same number of data points from the observed light curve allowing for multiple extraction. Then, the CCF and τ_cent are calculated for each realization in the same way to produce the cross-correlation centroid distribution (CCCD). The CCCDs calculated from 5000 realizations with the FR/RSS method are presented in Figure 10, and the derived τ_cent and its uncertainties are listed in Table 5. Thus, the results with σ_add = 0 mJy and 0.06 mJy agree well with each other. On average, the lag has been estimated as ${\tau }_{{\rm{cent}}}=-{2.02}_{-0.50}^{+0.55}$ days (1σ error), and the time lag is significantly non-zero, because the 99% confidence limit is τ < −0.68 days.

Table 5.  Time Lag of the X-Ray Flux Variation behind that of the B Band

Analysis	${\sigma }_{B,\mathrm{sys}}$ a	Time Lagb
	(mJy)	(day)
		1%	15.9%	50%	84.1%	99%
τcentc	0.00	−3.14	−2.56	−2.09	−1.54	−0.73
	0.06	−3.25	−2.49	−1.95	−1.41	−0.64
τJAVd	0.00	−4.42	−2.97	−1.90	−1.32	−0.66
	0.06	−4.07	−2.95	−2.12	−1.53	−0.90

Notes.

^aThe systematic error for the B-band photometries added to the flux errors in quadrature. ^bPercentiles. ^cThe time lag based on the CCF analysis. ^dThe time lag based on the JAVELIN software.

Download table as:  ASCIITypeset image

3.3.2. JAVELIN Analysis

The other method to estimate the lag is the JAVELIN software developed by Zu et al. (2011), which is widely employed not only in recent reverberation studies of the optical broad emission lines and the thermal dust emission of AGNs (e.g., Grier et al. 2012; Koshida et al. 2014; Peterson et al. 2014), but also in analysis of the lag between flux variations of their X-ray emission and the UV–optical continuum emission (McHardy et al. 2014; Shappee et al. 2014; Lira et al. 2015). It explicitly builds a model of a response light curve that is expressed as a convolution of a source light curve and a top-hat transfer function with a certain lag, and it fits the model to the data on flux variations in different bands, one being the source light curve and the others the response light curves. The source light curve is modeled as a stochastic process using a DRW model, and the posterior distributions of the time lag as well as other model parameters are estimated using the Bayesian Markov chain Monte Carlo method.

We use the B-band light curve as the source but the X-ray data as the response, again with τ > 0 meaning X-ray delays. The X-ray data of epochs 1–5 and the optical data of MJD = 56316–56511 were selected for the JAVELIN analysis, and σ_add = 0 or 0.06 mJy was added to the B-band photometric errors in quadrature, just as in the CCF analysis. The posterior distribution of the time lag calculated from 250,000 realizations of the JAVELIN software is presented in Figure 11, and the resultant τ_JAV and its uncertainties are also listed in Table 5. The two values of σ_add again gave very similar results. The time lags were estimated as ${\tau }_{{\rm{JAV}}}=-{1.90}_{-1.07}^{+0.58}$ days and ${\tau }_{{\rm{JAV}}}=-{2.12}_{-0.83}^{+0.59}$ days (1σ error) with σ_add = 0 and 0.06, respectively. Again, we can exclude the case of τ = 0, because the 99% confidence limit is τ < −0.78 days on average. In summary, the two methods have given consistent results.

Zoom In Zoom Out Reset image size

As shown in Table 3 and Figure 4(top), we found a significant flux variation in the X-ray reflection component as well, on a timescale of several months, although its variation in amplitude is much smaller than that of the primary X-rays. This indicates that the source region of the reflection component has an extent on a scale of ∼0.1 pc. Moreover, the flux variation of the reflection component lags slightly behind that of the hard PL component, on a scale of a week, or possibly longer because the X-ray sampling becomes very sparse after epoch 5.

According to the unified model of AGNs, the reflection component accompanied by the neutral Fe–Kα line is supposed to arise from a dust torus, but different origins such as outer accretion disks and the broad emission-line region have also been suggested (e.g., Awaki et al. 1991; Yaqoob & Padmanabhan 2004; Nandra 2006; Jiang et al. 2011; Gandhi et al. 2015; Minezaki & Matsushita 2015). Interestingly, the reverberation studies of NGC 3516 yielded a dust lag of ∼50–70 days (Koshida et al. 2014) and a lag of broad Balmer emission lines of ∼7–13 days (Peterson et al. 2004; Denney et al. 2010). These are comparable to the timescale of the variation of the X-ray reflection and its delay behind the hard PL component. In order to further examine the origin of the X-ray reflection component, a direct comparison of its variation with that of the dust emission and the broad Balmer emission lines would be necessary; this will be discussed in a forthcoming paper.

In the present study, we performed simultaneous X-ray and optical monitoring of the type 1.5 Seyfert galaxy NGC 3516 with Suzaku and Japanese ground-based telescopes—the Pirka, Kiso Schmidt, MITSuME, Nayuta, and Kanata telescopes. By applying spectral fitting to the X-ray data and differential image photometry to the B-band images, and quantitatively comparing the X-ray and optical flux variations, we have obtained the following results.

• 1.  
  During our observations, NGC 3516 was in an X-ray-faint state. It was faintest in epoch 7, when the 2–10 keV flux became ∼0.21 × 10⁻¹¹ erg s⁻¹ cm⁻², which is just 5% of the average flux recorded in 1997–2002 with RXTE. Even when brightest (epoch 4), the 2–10 keV flux was ∼2.70 × 10⁻¹¹ erg s⁻¹ cm⁻², which is only 70% of the average in 1997–2002. The flux varied on timescales longer than a few days, without any intraday changes.
• 2.  
  The 2–45 keV emission mainly consisted of two spectral components: a variable hard PL continuum with a photon index of ∼1.75 and a reflection component with a prominent narrow Fe–Kα emission line. Throughout the monitoring, the hard X-ray component retained almost the same spectral shape and exhibited a peak-to-peak flux change of ∼2.5 × 10⁻¹¹ erg s⁻¹ cm⁻², or about an order of magnitude.
• 3.  
  The B-band flux density varied by ∼4.0 mJy from peak to peak, which translates to a flux change of ∼2.7 × 10⁻¹¹ erg s⁻¹ cm⁻² in νF_ν units, after correcting for the optical extinction. The B-band flux varied on timescales longer than a few days, like that of the X-rays.
• 4.  
  X-ray and B-band flux correlation was significantly detected for the first time in NGC 3516. The flux changes of the hard X-ray component significantly preceded those in the B band by ${2.0}_{-0.6}^{+0.7}$ days (1σ error).

In previous monitoring of NGC 3516 conducted with RXTE and the Israeli WISE telescope in 1997–2002, no significant correlations were derived between X-ray and B-band flux variations, giving a cross-correlation coefficient (CCC) of <0.35 (Maoz et al. 2002). In contrast, we succeeded in detecting a significant correlation between them with a high CCC of >0.95. What is responsible for the clear difference? Maoz et al. (2002) argued that the lack of correlation might be due to changes in absorption that are independent of the primary variations in X-rays. However, the 2–10 keV flux varied during their monitoring by a factor of ∼4 from peak to peak. Such a large change would be hardly explained by the variations recorded so far in the column density of optically thick neutral absorbers (Turner et al. 2011). We hence suggest alternatively that the dominant X-ray-variable component changed between the two monitoring campaigns.

Noda et al. (2013a) discovered that the X-ray emission of NGC 3516 comprises at least two different primary continua with distinct spectral shape and flux variability; a flat continuum (Γ ∼ 1.1–1.7) and a steep one (Γ ∼ 2.3), which we hereafter call the hard and gradually varying primary component (HGPC) and the soft and rapidly varying primary component (SRPC), respectively. Noda et al. (2013a) showed that the luminosities of the HGPC and SRPC were comparable in a bright state in 2005, when the absorption-corrected 2–10 keV total flux was F_total ∼ 3.5 × 10⁻¹¹ erg s⁻¹ cm⁻². In contrast, the 2009 Suzaku observation caught NGC 3516 in a faint state with F_total ∼1.1 × 10⁻¹¹ erg s⁻¹ cm⁻², in which the variation was carried solely by a Γ = 1.7 PL, which can be identified with the HGPC. Thus, there is a certain threshold in between these two flux values, say, at F_th ∼ 3 × 10⁻¹¹ erg s⁻¹ cm⁻²; the HGPC and SRPC coexist above F_th, whereas only the HGPC remains below F_th. These two primary continua have also been identified in another Seyfert, NGC 3227 (Noda et al. 2014)¹⁹ : their luminosity-dependent behavior was found to be very similar to that observed from NGC 3516. Namely, the HGPC was always present, whereas the SRPC appeared when the source was brighter than a certain threshold, to then coexist with the other.

In the present Suzaku observations, the absorption-corrected 2–10 keV flux of NGC 3516 was in the range (0.2–2.7) ×10⁻¹¹ erg s⁻¹ cm⁻², i.e., always below the suggested F_th. Furthermore, the 2–45 keV spectrum of the main variable PL was flat with Γ ∼ 1.75, and the CCP distributions are smoothly connected to those in 2009 (Figure 2). Hence, we conclude that NGC 3516 was in the faint phase, with the HGPC dominating just as in the 2009 observation. On the other hand, during the 1997–2002 monitoring by Maoz et al. (2002), the averaged 2–10 keV flux was F_total ∼ 4.0 × 10⁻¹¹ erg s⁻¹ cm⁻², which is higher than F_th, indicating that NGC 3516 was mostly in the bright phase. If the SRPC has a significantly poorer correlation with the optical than the HGPC does, the correlation between the total X-ray flux and the optical should become worse when NGC 3516 is in the bright phase. Thus, the difference between the data of Maoz et al. (2002) and ours may be explained by presuming that the HGPC is well correlated with the optical while the SRPC is not.

To generalize the above conclusion, we collected CCCs from peak values of the CCF in previous simultaneous observations of various Seyferts. As a result, sources with high L_bol/L_Edd ratios were found to have relatively low CCCs: MCG–6-30-15 with L_bol/L_Edd ∼ 0.2 has CCC ≲ 0.1 (Lira et al. 2015), NGC 4051 with L_bol/L_Edd ∼ 0.15 has CCC ∼ 0.3 (Breedt et al. 2009), and MR 2251-178 with L_bol/L_Edd ∼ 0.2 has CCC ∼ 0.5 (Arévalo et al. 2008). In contrast, sources with low L_bol/L_Edd ratios possibly have high CCCs: NGC 6814 with L_bol/L_Edd ∼ 0.008 has CCC ∼ 0.9 (Troyer et al. 2016), and NGC 3516 in the present study has CCC ≳ 0.95 wherein L_bol/L_Edd < 0.01 is indicated by the value of 0.00612 (Vasudevan & Fabian 2009), obtained in 2001 when NGC 3516 showed F_tot close to the second highest among ours. Therefore, we speculate that X-ray–optical correlation in Seyfert galaxies becomes commonly worse when the Eddington ratio increases, because the SRPC starts to dominate.

According to Figure 12(b) and related descriptions in Noda et al. (2014), the SRPC may originate from patchy coronae heated by a magnetic process on the surface of a disk (e.g., Reynolds & Nowak 2003). If the patchy coronae cover a small fraction of the disk, the SRPC generated in the coronae can be related to just small areas of the disk, making X-ray and optical variations independent. Another possible origin of the SRPC is a faulty jet at the accretion axis of the SMBH (e.g., Ghisellini et al. 2004). If the region of jet emission is located far away from the accretion disk as suggested by Noda et al. (2014), the connection between the SRPC and the disk optical emission may become weaker, making their correlation poorer. In any case, the emergence and disappearance of these SRPC-generating regions are likely to be more localized effects than those responsible for the transitions from hard to soft state seen in black-hole binaries, because the changes in AGN state considered here take place on timescales of weeks to years, which are much shorter than would be expected (e.g., several tens of thousands of years) if the typical state-transition timescales of a few days are scaled to the mass ratios.

Zoom In Zoom Out Reset image size

4.3.1. Application of the Standard X-Ray Reprocessing Model

The measured optical lag of ${2.0}_{-0.6}^{+0.7}$ days strongly supports the X-ray reprocessing model (e.g., Krolik et al. 1991; Cackett et al. 2007); optical brightening occurred, at least in the present case, through irradiation by the increased hard X-ray intensity. This process is energetically feasible, because the increase in 2–10 keV X-ray flux by 2.5 × 10⁻¹¹ erg s⁻¹ cm⁻² is comparable to the increase in B-band flux by ∼2.7 × 10⁻¹¹ erg s⁻¹ cm⁻². The optical time lag in NGC 3516 obtained here is similar to those of a few days observed in other Seyfert galaxies, including NGC 5548 (e.g., Suganuma et al. 2006; McHardy et al. 2014; Edelson et al. 2015), NGC 3783 (Arévalo et al. 2009), NGC 4051 (Breedt et al. 2010), and NGC 6418 (Troyer et al. 2016). Therefore, the suggested mechanism, i.e., X-ray reprocessing, can be operating commonly among these Seyfert galaxies.

Let us, then, investigate whether or not the derived optical lag can be explained by the most commonly adopted X-ray reprocessing model (e.g., Cackett et al. 2007), which assumes that an optically thick, geometrically thin, standard accretion disk continues down to the innermost last stable circular orbit (ISCO), located at $3{R}_{{\rm{s}}}$ where ${R}_{{\rm{s}}}=2{{GM}}_{{\rm{BH}}}/{c}^{2}$ is the Schwarzschild radius and M_BH the black-hole mass, and is illuminated by a lamppost-type X-ray source proximate to the black hole. At radii R ≫ R_s, the radial temperature profile of the disk is given as (Shakura & Sunyaev 1973)

$\begin{eqnarray}&&T(R)={\left(\displaystyle \frac{3{{GM}}_{{\rm{BH}}}\dot{M}}{8\pi \sigma }\right)}^{\tfrac{1}{4}}{R}^{-\tfrac{3}{4}},\end{eqnarray} \tag{ 1 }$

where $\dot{M}$ is the mass accretion rate, G is the gravitational constant, and σ is the Stefan–Boltzmann constant.

After a delay time τ, the X-ray flux variation that arose close to the black hole will propagate to a radius $R=c\tau$, and will increase the emissivity there. Let us define a wavelength $\lambda ={hc}/{kT}(c\tau )X$ in such a way that the continuum emission at this radius is peaked at λ, where X is a numerical factor of order unity. Although Wien's displacement law for blackbody radiation in νB_ν gives X = 3.92, we here set X ∼ 3.2, according to the detailed calculations of the disk response to continuum emission by Cackett et al. (2007) and Collier et al. (1999). Substituting $c\tau$ for R, and expressing T(cτ) with λ, Equation (1) can be rewritten as

$\begin{eqnarray}&&T(R)=\displaystyle \frac{{hc}}{k\lambda X}{\left(\displaystyle \frac{R}{c\tau }\right)}^{-\tfrac{3}{4}}.\end{eqnarray} \tag{ 2 }$

By eliminating $T(R)\,{R}^{3/4}$ from Equations (1) and (2), we obtain $\dot{M}\propto {M}^{-1}{\tau }^{3}{\lambda }^{-4}$. Further expressing $\dot{M}$ as the bolometric luminosity of the accretion disk, ${L}_{{\rm{bol}}}\sim 0.1\,\dot{M}{c}^{2}$ (assuming that the disk extends to the ISCO), and normalizing it to the Eddington luminosity (∝ M_BH), the Eddington ratio of the source η is estimated as

$\begin{eqnarray}\eta & \equiv & {L}_{{\rm{bol}}}/{L}_{{\rm{Edd}}}\sim 3.8\ {\left(\displaystyle \frac{{M}_{{\rm{BH}}}}{3.2\times {10}^{7}{M}_{\odot }}\right)}^{-2}\\ & & \times {\left(\displaystyle \frac{\tau }{2.0{\rm{day}}}\right)}^{3}{\left(\displaystyle \frac{\lambda }{0.44\mu {\rm{m}}}\right)}^{-4}.\end{eqnarray} \tag{ 3 }$

Equation (3), in fact, implies a serious problem. First of all, it indicates an unrealistic "super-Eddington" condition. Even putting aside this issue, it means a discrepancy of a factor of ∼500 from the likely value of η = 0.006–0.01 (Section 4.2), which characterizes the nucleus of NGC 3516 during the present observations. This difficulty can be stated more directly: by summing up the blackbody contributions from various disk annuli, the continuum spectrum from the accretion disk under consideration can be calculated as

$\begin{eqnarray}{f}_{\nu } & \sim & 550\ {\rm{mJy}}\ {\left(\displaystyle \frac{\tau }{2.0{\rm{days}}}\right)}^{2}{\left(\displaystyle \frac{D}{41.3{\rm{Mpc}}}\right)}^{-2}\\ & & \times {\left(\displaystyle \frac{\lambda }{0.44\mu {\rm{m}}}\right)}^{-3}\cos i,\end{eqnarray} \tag{ 4 }$

where i is the disk inclination (Collier et al. 1999; Cackett et al. 2007). Thus, the measured delay predicts very bright optical emission, arising from the inner part of the accretion disk. In other words, the SMBH would have to be accreting at a very high rate in order for its accretion disk to emit the B-band light from such a large distance as ∼2 lt-day. Equation (4) means a contradiction by two orders of magnitude with the present optical data, in which the peak-to-peak flux variation in the B band was ∼4.0 mJy, and the minimum B-band flux during the observation would be less than that. Conversely, if we started from the assumption of η = 0.006–0.01, the B-band lag predicted by Equation (3) would become ∼8 times smaller than was measured. It is thus extremely difficult to reconcile the measured clear optical delay of ∼2 days with the observed optical faintness, as long as we assume a standard disk extending down to the ISCO and an illuminating X-ray source close to the SMBH.

4.3.2. Consideration of the X-Ray Irradiation

In the X-ray reprocessing model, we should consider not only the viscous heating but also the heating by the X-ray irradiation. The former, which underlies Equation (1), is written as

$\begin{eqnarray}&&{D}_{{\rm{vis}}}=\displaystyle \frac{3{{GM}}_{{\rm{BH}}}\dot{M}}{8\pi {R}^{3}},\end{eqnarray} \tag{ 5 }$

whereas the latter is given as

$\begin{eqnarray}&&{D}_{{\rm{irr}}}=\displaystyle \frac{(1-A){L}_{{\rm{X}}}}{4\pi {R}_{{\rm{X}}}^{2}}\cos \theta ,\end{eqnarray} \tag{ 6 }$

where L_X is the X-ray luminosity, A is the disk albedo, R_X is the distance from the disk to the X-ray illuminator, and θ is the angle of incidence of X-rays onto the disk. Collier et al. (1999) and Cackett et al. (2007) assumed that the X-ray emitter is located on the rotating axis of the SMBH like a "lamppost" (Figure 12(a)) and its height H_X from the SMBH is much smaller than R (H_X ≪ R). Under these assumptions, we obtain R_X ∼ R and $\cos \theta ={H}_{{\rm{X}}}/\sqrt{{R}^{2}+{H}_{{\rm{X}}}^{2}}\sim {H}_{{\rm{X}}}/R$. Then, D_vis and D_irr have just the same R-profile, and the total heating per unit face area of the disk, ${D}_{{\rm{tot}}}={D}_{{\rm{vis}}}+{D}_{{\rm{irr}}}$, is expressed as

$\begin{eqnarray}&&{D}_{{\rm{tot}}}=\displaystyle \frac{3{{GM}}_{{\rm{BH}}}\dot{M}}{8\pi {R}^{3}}+\displaystyle \frac{(1-A){L}_{{\rm{X}}}{H}_{{\rm{X}}}}{4\pi {R}^{3}}.\end{eqnarray} \tag{ 7 }$

As a result, the radial temperature profile of the disk is now given by

$\begin{eqnarray}&&T(R)={\left[\displaystyle \frac{3{{GM}}_{{\rm{BH}}}}{8\pi \sigma }\left\{\dot{M}+\displaystyle \frac{2(1-A){L}_{{\rm{X}}}{H}_{{\rm{X}}}}{3{{GM}}_{{\rm{BH}}}}\right\}\right]}^{\tfrac{1}{4}}{R}^{-\tfrac{3}{4}}.\end{eqnarray} \tag{ 8 }$

This exhibits the same R-dependence as the standard viscous disk, namely, T ∝ R^−3/4. Furthermore, because of the proportionality ${L}_{{\rm{X}}}\propto \dot{M}$, we should use $\dot{M}+2(1-A){L}_{{\rm{X}}}{H}_{{\rm{X}}}/3{{GM}}_{{\rm{BH}}}$ in place of $\dot{M}$ when converting this expression into that of η according to Equation (3). Consequently, the Eddington ratio η and the B-band flux density can be given by the same Equations (3) and (4), respectively.²⁰ Thus, we encounter just the same problem as before, even when the X-ray irradiation is taken into account.

According to, e.g., Cackett et al. (2007), McHardy et al. (2014), and Edelson et al. (2015), typical Seyfert galaxies including NGC 5548, NGC 4051, and Mrk 335 yielded values of τ that are a few times larger than those expected from their estimated η. According to Morgan et al. (2010), the accretion disks of gravitationally lensed quasars are larger, by a factor of ∼4, than those expected from η, just like in the Seyfert galaxies. These subtle but systematic contradictions between the disk size and η have also been well known in the X-ray reprocessing scenario (e.g., Collier et al. 1999; Cackett et al. 2007). Recently, Troyer et al. (2016) reported that NGC 6814 exhibited a value of τ that is ∼40 times larger than that predicted from the estimated η. In the present case of NGC 3516, the difference in the measured versus predicted optical delay amounts to a factor of ∼10, which is considerably larger than those of NGC 5548, NGC 4051, and Mrk 335. These systematic discrepancies imply that the lamppost-type geometry of X-ray irradiation shown in Figure 12(a) has a common problem, which is considered to become more prominent toward lower luminosities, considering that NGC 3516 and NGC 6814 have both relatively low values of η. Thus, we need to revise the geometrical assumptions in order to make τ and η consistent in Seyfert galaxies with low Eddington ratio.

4.3.3. Possible Geometry of Accretion Flows in Seyferts with Low Eddington Ratio

The present coordinated X-ray and optical observations of NGC 3516 found the object in a very faint state at both wavelengths. A tightly correlated change in intensity lasting for ∼60 days was detected, wherein the optical B-band variation showed a clear delay by ∼2 days behind that of the X-ray continuum, which is described by a PL model with Γ ∼ 1.7. This optical lag indicates the effect of X-ray reprocessing, but the delay cannot be reconciled with the B-band faintness as long as we consider the standard reprocessing scenario, which invokes a standard accretion disk irradiated by a lamppost-type X-ray source. Instead, the observed results can be consistently explained by assuming that the disk is truncated at ∼500 R_s and turns into a RIAF where the illuminating hard X-rays are produced via thermal Comptonization.