A complex dust morphology in the high-luminosity AGN Mrk 876

Recent models for the inner structure of active galactic nuclei (AGN) advocate the presence of a radiatively accelerated, dusty outflow launched from the outer regions of the accretion disk. Here we present the first near-infrared (near-IR) variable (rms) spectrum for the high-luminosity, nearby AGN Mrk 876. We find that it tracks the accretion disk spectrum out to longer wavelengths than the mean spectrum due to a reduced dust emission. The implied outer accretion disk radius is consistent with the infrared results predicted by a contemporaneous optical accretion disk reverberation mapping campaign and much larger than the self-gravity radius. The reduced flux variability of the hot dust could be either due to the presence of a secondary, constant dust component in the mean spectrum or introduced by the destructive superposition of the dust and accretion disk variability signals or some combination of both. Assuming thermal equilibrium for optically thin dust, we derive the luminosity-based dust radius for different grain properties using our measurement of the temperature. We find that in all cases considered the values are significantly larger than the dust response time measured by IR photometric monitoring campaigns, with the least discrepancy present relative to the result for a wavelength-independent dust emissivity law, i.e. a blackbody, which is appropriate for large grain sizes. This result can be well explained by assuming a flared, disk-like structure for the hot dust.


INTRODUCTION
Most active galactic nuclei (AGN) display prominent infrared (IR) emission in their spectral energy distributions (SEDs), which can be attributed to thermal dust radiation.This central dusty structure is commonly assumed to be optically thick and have a toroidal geometry Corresponding author: Hermine Landt hermine.landt@durham.ac.uk aligned with the plane of the accretion disk.The latter requirement was mainly derived from the need of an equatorial obscurer in AGN in order to account for observations of broad emission lines in the polarized, scattered light of many type 2 AGN and the relative numbers of type 1 and type 2 AGN in complete samples (see reviews by Lawrence 1987;Antonucci 1993;Netzer 2015;Lyu & Rieke 2022a).A warped, dusty disk could be a viable alternative to the dusty torus, as first proposed by Phinney (1989).The extended dusty torus emits over a large range of IR wavelengths, with the hottest dust observed in the near-IR believed to be located closest to the central supermassive black hole, whereas warm dust (of a few 100 K) is observed in the mid-IR and expected to be located on average further out (Weedman et al. 2005;Buchanan et al. 2006;Landt et al. 2010;Lyu & Rieke 2018;Brown et al. 2019).Often, cold dust (of a few 10 K) emitting at far-IR wavelengths is also observed in AGN, although this component is most likely heated by young stars in the host galaxy rather than by the ultraviolet (UV) and optical emission of the central accretion disk (Kirkpatrick et al. 2015).
AGN are ideal laboratories to investigate the chemical composition and grain properties of astrophysical dust.Their UV/optical luminosities are usually high enough to heat the circumnuclear dust to sublimation temperatures.If these highest values can be observed, they can in principle constrain the chemistry since different species condense out of the gas phase in different environmental conditions.Dust temperatures were measured with simultaneous photometry at several near-IR wavelengths in a handful of sources (Clavel et al. 1989;Glass 2004;Schnülle et al. 2013Schnülle et al. , 2015)), but has only come of age with the availability of efficient near-IR cross-dispersed spectrographs.Landt et al. (2011b) and Landt et al. (2014) derived dust temperatures from such spectroscopy for the largest sample of type 1 AGN so far (∼ 30 sources).Their measurements yielded a very narrow temperature distribution with an average value of T ∼ 1400 K.This result indicated that, if the hot dust is indeed at sublimation, it is composed of only silicate dust grains and so an oxygen-rich environment from which the dust formed.However, if carbon is present, then the dust is not heated to close to sublimation, since carbonaceous dust, e.g.graphite, can survive up to T ∼ 2000 K (Salpeter 1977).
The range in temperatures within the dusty torus roughly translates to a range in radius, whereby mainly the mid-IR emitting dust can be spatially resolved with current instruments.High-angular-resolution mid-IR imaging (Ramos Almeida et al. 2011) and mid-IR interferometry (e.g.Tristram et al. 2009;Burtscher et al. 2013) of a dozen nearby and bright AGN delivered useful upper limits on the extent of the dusty torus (of up to a few parsec).Such spatially resolved observations have also revealed a significant warm dust component perpendicular to the plane of the accretion disk, referred to as 'polar dust' (Hönig et al. 2013;Tristram et al. 2014;Isbell et al. 2022;Lyu & Rieke 2022b;Gámez Rosas et al. 2022), which was also evident from SED studies (Landt et al. 2010;Isbell et al. 2021).GRAVITY, the near-IR interferometric instrument at the Very Large Telescope (VLT), has started to resolve the innermost and hottest part of the central dusty structure in some nearby, luminous sources (Gravity Collaboration et al. 2020) and further progress is expected from an upgrade in sensitivity to GRAVITY+.But for most AGN knowledge about the inner dust radius is most efficiently obtained through dust reverberation, which is a technique that measures the time response of the dust to the variable, irradiating accretion disk flux.For ∼ 60 AGN, the hot dust radius was determined by photometric campaigns, often coordinated at optical and near-IR wavelengths (e.g.Clavel et al. 1989;Nelson 1996;Oknyanskij & Horne 2001;Glass 2004;Suganuma et al. 2006;Schnülle et al. 2013;Koshida et al. 2014;Schnülle et al. 2015;Vazquez et al. 2015;Minezaki et al. 2019;Lyu et al. 2019).In general, observed dust response times follow a luminosity-radius relationship with a slope similar to that for the broad emission line region (BLR), indicating a roughly constant hot dust flux and so narrow hot dust temperature distribution.However, dust radii measured via reverberation are often smaller (by factors of a few) than dust radii measured by interferometry or estimated from the SED assuming thermal equilibrium (Kishimoto et al. 2007;Nenkova et al. 2008;Landt et al. 2014).A possible interpretation for this finding is that the dust has a bowl-shaped geometry caused by the anisotropy of the accretion disk emission irradiating it (Kawaguchi & Mori 2010, 2011).In such a geometry, the hot dust located in the plane of the accretion disk is placed further in than the bulk of the dust and so is expected to dominate the reverberation signal.
Progress in technology and flexible scheduling have now made spectroscopic near-IR monitoring campaigns feasible.Landt et al. (2019) presented the first such program.Using a near-IR cross-dispersed spectrograph with a relatively wide wavelength coverage that extends partially into the optical they showed that such a campaign allows in low-redshift AGN: (i) the monitoring of a large portion of the hot dust SED that readily gives the dust temperature and its evolution; (ii) the separation of the accretion disk emission from that of the hot dust, which is a considerable source of error in photometric campaigns; (iii) the simultaneous determination of luminosity-based and response-weighted dust radii; (iv) the construction of the variable (rms) near-IR spectrum; and (v) the study of the variability of emission lines formed in the BLR and the coronal line region.Their study of the hot dust in NGC 5548 found that a single component dominated both the mean emission and variations, with the dust reponse time and the luminosity-based dust radius being consistent with each other only if a blackbody emissivity was assumed.This result constrained the dust grain size to a few µm.The temperature and its variability indicated carbonaceous dust well below the sublimation threshold undergoing a heating and cooling process in response to the variable UV/optical accretion disk flux irradiating it.Most importantly, the dust reverberation signal showed tentative evidence for a second hot dust component most likely associated with the accretion disk.The existence of such dust is a prerequisite for the recent models of the AGN structure proposed by Czerny et al. (2017) and Baskin & Laor (2018), which explain both the BLR and dusty torus as part of the same outflow launched from the outer regions of the accretion disk by radiation pressure on dust, preferentially on carbonaceous dust since it has a higher opacity than silicate dust.
Here we present results from a near-IR spectroscopic monitoring campaign on Mrk 876 conducted between 2016 May and 2017 Jul.A contemporaneous optical photometric monitoring campaign in several bands was conducted between 2016 Mar and 2019 May and Miller et al. (2022) have recently presented the accretion disk reverberation results from these data.The structure of our paper is as follows.In Sections 2, we discuss the science target and give the details of the near-IR observations and data reduction in Section 3. In Sections 4 and 5, we derive the luminosity-based dust radius and the variable (rms) spectrum, respectively.We discuss our main results in Section 6, where we seek a relation between the central dust structures observed in AGN and protoplanetary disks around young stars.Finally, in Section 7, we present a short summary and our conclusions.

THE SCIENCE TARGET
Mrk 876 (PG 1613+658) is one of the intrinsically most luminous AGN in the nearby Universe.At a redshift of z = 0.129 it has an average V -band luminosity of ∼ 5 × 10 44 erg s −1 (Bentz et al. 2013).This places it at the undersampled top end of the relationship between the hot dust radius and optical continuum luminosity presented by Koshida et al. (2014) and Minezaki et al. (2019).The latter study performed for Mrk 876 a co-ordinated optical and near-IR photometric campaign during the years 2003-2007 and measured a dust response time of τ = 320 ± 34, 327 +42 −36 and 327 +41 −35 lightdays, assuming different values for the accretion disk power-law spectral slope when decomposing the flux in the K band.A similar dust response time was obtained also by Lyu et al. (2019), who combined optical photometric light-curves from ground-based transient surveys with the WISE monitoring data for the years 2010-2018.They estimated a dust response time of τ = 372 ± 63 and 390 ± 37 light-days for the W 1 and W 2 bands, re-spectively.Finally, Afanasiev et al. (2019) performed spectropolarimetric observations of Mrk 876 in 2015 and constrained the distance to the equatorial scattering region to 254 ± 39 light-days, which they assumed was the inner radius of the dusty torus.
Mrk 876 has a relatively high black hole mass, well-determined by optical reverberation campaigns of M BH = (2.2 ± 1.0) × 10 8 M (Peterson et al. 2004;Bentz & Katz 2015), transformed from the measured virial product using a scaling factor of f = 4.3 (Grier et al. 2013).For this black hole mass, the corresponding gravitational radius is r g = 3.2 × 10 13 cm = 0.013 lightdays, where r g = GM BH /c 2 , with G the gravitational constant and c the speed of light.The corresponding Eddington luminosity is L edd = 2.7 × 10 46 erg s −1 .The optical emission-line spectrum of Mrk 876 is of the inflected type, i.e. the broad line profiles have clearly discernible broad-and narrow-line components.The hydrogen BLR radius has been recently measured by an optical spectroscopic reverberation mapping campaign during the years 2016-2021 to be in the range of ∼ 40 − 50 lightdays, based on the response time of the optical Balmer line Hβ (Bao et al. 2022).
Mrk 876 (J2000 sky coordinates R.A. 16 h 13 m 57.2 s , Decl.+65 • 43 10 ) is observable only from the northern hemisphere.Its low redshift and high intrinsic luminosity make it sufficiently bright in the near-IR (2MASS J = 13.2 mag, K s = 11.4 mag; Skrutskie et al. 2006) so that a high-quality spectrum can be achieved in a relatively short exposure time.We adopt the cosmological parameters H 0 = 70 km s −1 Mpc −1 , Ω M = 0.3, and Ω Λ = 0.7, which give a luminosity distance to Mrk 876 of 605.2 Mpc and an angular scale at the source of 2302 pc per arcsec.

The near-IR spectroscopy
We observed Mrk 876 between 2016 May and 2017 Jul at the Gemini North 8 m observatory on Maunakea, Hawaii, in queue mode (Program ID: GN-2016A-FT-27, GN-2017A-Q-41) using the Gemini Near-Infrared Spectrograph (GNIRS; Elias et al. 2006).We obtained roughly one spectrum per month, except for the gap period of about five months (Sep to Jan) when Mrk 876 is not observable, resulting in a total of 10 near-IR spectra.Table 1 lists the journal of observations.The sky was clear and the seeing was very good (∼ 0.4 − 0.6 ) for all nights except 2016 Jun 16, 2017 Feb 24 and Jul 5, when some clouds were present and the seeing was variable.We used the cross-dispersed mode with the short camera at the lowest spectral resolution (31.7 l mm −1 grating), thus covering the entire wavelength range of 0.85 − 2.5 µm without inter-order contamination.We chose a slit of 0.675 × 7 , which we oriented at the parallactic angle.This set-up gives an average spectral resolution of full width at half-maximum (FWHM) ∼ 400 km s −1 .The on-source exposure time of 8 × 120 s ensured that we obtained spectra with a high signalto-noise ratio (S/N ) in order to reliably measure the emission-line profiles.Since the source is not extended in the near-IR, we nodded along the slit in the usual ABBA pattern.The observations were taken on average at an airmass of sec z ∼ 1.5, which is close to the minimum value achievable for this source from Hawaii.
After or before the science target, we observed the nearby (in position and airmass) F0 V star HIP 79087 (HD 145710) that has accurate near-IR magnitudes.We used this standard star to correct our science spectrum for telluric absorption and for flux calibration.For the flux calibration we assumed that its continuum emission can be approximated by a blackbody with an effective temperature of T eff = 6769 K (Soubiran et al. 2016).Flats and arcs were taken after the science target.
We reduced the data using the Gemini/IRAF package with GNIRS specific tools (Cooke & Rodgers 2005).The data reduction steps included preparation of calibration and science frames, processing and extraction of spectra from science frames, wavelength calibration of spectra, telluric correction and flux-calibration of spectra, and merging of the different orders into a single, continuous spectrum.The spectral extraction width was adjusted interactively for the telluric standard star and the science source to include all the flux in the spectral trace.A local averaged background flux was fitted and subtracted from the total source flux.The final spectrum was corrected for Galactic extinction using the IRAF task onedspec.dereddenwith an input value of A V = 0.005, which we derived from the Galactic hydrogen column densities published by Dickey & Lockman (1990).In Fig. 1, we show the spectrum from 2016 May 25 as a representative example.

The complementary photometry
Mrk 876 was monitored with the 1 m robotic telescope network of the Las Cumbres Observatory (LCOGT; Brown et al. 2013) as part of the 2014 AGN Key Project almost daily in four bands (g, r, i and z s ) between 2016 Feb 12 and 2020 Oct 10, but with the z s band switched off in the period of 2017 Apr 4 to 2018 Jun 17.We made use here of the light-curves in the g and z s filters, which have a wavelength width of 1500 Å and 1040 Å around their central wavelength of 4770 Å and 8700 Å, respectively.The LCOGT Sinistro cameras have a FOV of 26.5 × 26.5 and a plate scale of 0.389 per pixel.Gemini GNIRS near-IR spectrum from 2016 May 25, shown as observed flux versus rest-frame wavelength.Emission lines listed in Table 4 of Landt et al. (2008) are marked by dotted lines and labeled; black: permitted transitions, green: permitted Fe II multiplets (not labeled), red: forbidden transitions and cyan: forbidden transitions of iron (those of [Fe II] not labeled).
The frames were first processed by LCOGT's BANZAI pipeline (McCully et al. 2018) in the usual way (bias and dark subtraction, flat-fielding correction and cosmic ray rejection) and were subsequently analysed with the custom-made pipeline described in detail by Hernández Santisteban et al. (2020).In short, after performing aperture photometry with a diameter of 7 and subtracting a background model, stable light-curves were produced by constructing a curve of growth using the standard stars on each individual frame and measuring the correction factors required to bring all different light-curves to a common flux level.A colour-correction and correction for atmospheric extinction were applied before the photometric calibration.Finally, an image zero-point calibration was performed at each epoch based on comparison stars in the field.In Fig. 2 (top two panels), we display the g and z s light curves, with the latter overlapping in wavelength with the GNIRS near-IR spectrum.

The absolute spectral flux scale
In order to derive light-curves, construct meaningful mean and variable (rms) spectra and estimate luminosities, we must achieve an accurate absolute flux calibration of the spectra.Similar to the approach chosen by optical spectroscopic reverberation campaigns to achieve this, we base our absolute spectral flux scale on a strong narrow emission line from a forbidden transition.Such emission lines are expected to remain constant during the campaign, since they are produced in gas that is located at distances of pc-to-kpc scales from the ionising source and has a number density low enough for recombination timescales to be large.We have chosen the [S III] λ9531 line since it is the strongest narrow forbidden emission line in the near-IR.It is blended with the Pa broad emission line, but can be easily separated from it since this hydrogen line is usually weak.Imaging studies of the [S III] λ9531 line have been performed in a few nearby AGN using near-IR integral-field spectrographs and the gas forming this line was found to be extended on moderate scales of up to a few 100 pc (e.g.Storchi-Bergmann et al. 2009;Fischer et al. 2015Fischer et al. , 2017)).Our slit size of 0.675 centered on the nucleus corresponds to a radial extent of ∼ 800 pc, which is sufficiently large to ensure that all the [S III] line flux is enclosed in the spectral extraction aperture.
We have used two routines to determine the photometric correction factors for our observations, the Prep-Spec routine developed by Keith Horne (last described in Horne et al. 2020) and the mapspec routine of Fausnaugh (2017).In short, PrepSpec models both the emission lines and the total continuum and subsequently matches the profiles of selected narrow emission lines in the spectra in order to derive the time-dependent scaling factor.It assumes that the majority of the spectra are photometric and so holds the median photometric correction factor at one.The mapspec routine constructs a reference spectrum by averaging the highest quality spectra, e.g., Gemini GNIRS near-IR spectrum from 2016 May 25 shown as luminosity versus rest-frame wavelength (black spectrum).We decomposed the continuum into two components, namely, an accretion disk spectrum that approximates the wavelength range ∼ < 1 µm (black dotted curve) and still dominates at 1 µm (vertical dashed line) and hot dust emission (red spectrum).We fitted the hot dust continuum with a blackbody spectrum (red dotted curve) and modified blackbody spectra for carbon and silicate dust (green and blue dotted curves, respectively) with resulting best-fit temperatures as listed in Table 2.
those spectra that have the highest S/N and/or were observed under the best weather conditions.Then, the profile of the chosen narrow emission line is matched between the reference spectrum and the individual spectra.The spectral scaling factors resulting from the two routines are listed in Table 1.With the exception of the last observational epoch, the results from PrepSpec and mapspec are consistent within the errors, but there is a trend for the former correction factors to be lower than the latter.In the following, we use the PrepSpec correction factors, which span a range of ∼ 1 − 40%, with an average value of ∼ 10%.In Fig. 2 (bottom two panels), we show the GNIRS near-IR spectral continuum flux in the observed wavelength range of the LCOGT z s band filter centred at 8700 Å, both original and corrected.If we subtract from the z s band light-curve a constant host galaxy flux contribution of ∼ 20% of the total flux, it matches very well all four of the corrected near-IR spectral flux measurements in the time of overlap.For the later epochs, the g band light-curve shows a downward trend in flux, which is also reproduced.In order to be able to derive luminosity-based dust radii, we need to measure the dust temperature and estimate the UV/optical accretion disk luminosity that heats the dust.The relatively large wavelength range of the cross-dispersed near-IR spectra covers in Mrk 876 simultaneously about half the hot dust SED and a considerable part of the accretion disk spectrum, which is expected to dominate the total continuum flux up to ∼ 1 µm (Landt et al. 2011a,b).Therefore, we decomposed the spectral continuum into two components following the approach described in Landt et al. (2019) and exemplified in Fig. 3.We note that, since we used a relatively small spectral aperture, the contribution of host galaxy light to the total observed continuum flux is negligible in this luminous AGN.In short, we have first approximated the rest-frame wavelength range of ∼ < 1 µm with the spectrum of a standard accretion disk, which we have subsequently subtracted from the total spectrum.For the calculation of the accretion disk spectrum we adopted the black hole mass given in Section 2 and the accretion rate was obtained directly from the a scaling to the data.Furthermore, we assumed that the disk is relatively large and extends out to r out = 10 4 r g .We then fitted the resultant hot dust spectrum at wavelengths > 1 µm with a blackbody, representing emis-sion by large dust grains, and also with two blackbodies modified by a power-law of the form Q λ (a) ∝ λ β , approximating with β = −1 and β = −2 the emissivity of sub-micron silicate and carbon dust grains, respectively (Landt et al. 2019, see their Fig. 8).We note in Fig. 3 the excellent agreement between the data and the standard accretion disk spectrum that we assumed here, a result generally found for AGN with negligible contribution from host galaxy light (e.g.Koratkar & Blaes 1999;Landt et al. 2011b).Table 2 lists the physical parameters derived from the spectral decomposition.The resultant average temperatures and the error on the mean are T = 1306 ± 6, 1134±5 and 1033±4 K for an emissivity law with β = 0, −1 and −2, respectively.
As in Landt et al. (2019), we then calculated luminosity-based dust radii, R d,lum , from the best-fit dust temperatures assuming radiative equilibrium between the luminosity of the irradiating source and the dust: where σ is the Stefan-Boltzmann constant and Q em is the Planck-averaged value of Q λ (a).We have approximated L uv with the accretion disc luminosity at each epoch, as given in Table 2, and have used for the Planck-averaged emission efficiencies in the case of β = −1 a value of Q em = 0.0210 appropriate for silicates of T = 1259 K and a = 0.1 µm (Laor & Draine 1993) and in the case of β = −2 a value of Q em = 0.0875 appropriate for graphite of T = 1000 K and a = 0.1 µm (Draine 2016).The average luminositybased rest-frame dust radii and the error on the mean are R d,lum = 469 ± 30, 4296±273 and 2694±172 lightdays in the case of a blackbody, and small-grain silicate and carbon dust, respectively.

THE VARIABLE NEAR-IR SPECTRUM
We calculated the mean and variable (rms) near-IR spectra following Peterson et al. (2004).Fig. 4 (top panel) shows them in comparison.In order to be able to compare the chromatic shape of the two spectra, we rescaled the rms spectrum by a factor of ∼ 7 in order to match the mean optical flux at 0.9 µm.Assuming a simple flux variability, scaling with the total flux would lead to comparable shapes, but we observe this similarity only in the spectral part dominated by the accretion disk.Instead two spectral continuum features become apparent at wavelengths beyond ∼ 1 µm.First, the rms spectrum tracks the accretion disk spectrum out to longer wavelengths (∼ 1.2 µm) than the mean spectrum (∼ 1 µm).This is an important finding, since we do not know the full extent of accretion disks, mainly because their infrared spectrum is overwhelmed by the emission from the hot dust.Infrared polarimetry can reveal the disk near-IR continuum, but only for the brightest AGN (Kishimoto et al. 2008).In order for the accretion disk spectrum to dominate the emission at ∼ 1 µm, as observed in the mean spectrum and established by the near-IR radius-luminosity relationship (Landt et al. 2011a), the inferred disk outer radius is r out ∼ 1500 r g ∼ 20 light-days.Then, observing the optical variable component in the rms spectrum out to ∼ 1.2 µm and assuming it is the accretion disk implies a value of r out ∼ 1900 r g ∼ 25 light-days.This radius is much larger than the self-gravity radius, i.e. the radius where accretion disks are generally assumed to become unstable and fragment, which is ∼ 12 light-days (Lobban & King 2022).
Secondly, the near-IR flux variability relative to the mean continuum flux at ≥ 1 µm appears significantly lower than in the visible domain, which is dominated by the accretion disk only.Since the near-IR shows clearly the additional hot dust flux, this behaviour could be explained by a systematically reduced flux variability of the hot dust with respect to the accretion disk.This would explain why the accretion disk spectral slope is revealed out to longer wavelengths, a situation similar to that seen in the mean flux spectrum of the so-called "hot-dust-poor" quasars reported by Hao et al. (2010).A spectral decomposition of the rms and mean spectra, as performed in Section 4, gives a flux ratio reduced by a factor of ∼ 3 in the former relative to the latter.The ratio of the rms to the mean spectrum measures the fractional amplitude of the variations, i.e. the variable flux as a percentage of the mean flux, in dependence of wavelength.If the shape of the variable spectrum is identical to that of the mean spectrum, it means that the fractional variability is independent of wavelength.A reduction of the fractional variability at a given wavelength can be introduced by the presence of a constant component, which would increase the flux in the mean spectrum.Our observed variable near-IR spectrum for Mrk 876 shows a reduced fractional variability for the blackbody component only (Fig. 5, bottom panel).We visualize this difference spectrum for wavelengths ∼ > 1 µm, which can be well-fitted by a blackbody of temperature T ∼ 1300 K, in Fig. 5 (top panel).Since the host galaxy of Mrk 876 is unusually luminous (e.g.Bentz et al. 2009), we have also considered the case that some of the spectral difference between the mean and rms spectrum is induced by this constant component.In Fig. 5 (top panel), we have included the template of an elliptical galaxy.It is clear that a significant contribution from this component is unlikely given that the flux difference between the mean and rms spectrum increases with wavelength, whereas the host galaxy spectral flux shows the opposite behaviour.
Therefore, the hot dust seems to be composed of at least two components with very different variability timescales or behaviours.The secondary dust component present in the mean spectrum appears constant in our campaign either because it is truly non-variable or because it varies on timescales not probed by our monitoring campaign.If such a constant dust component is indeed present, then the variable dust component does not seem to dominate the total dust emission and constitutes only ∼ 30% of it.A spectral decomposition of the variable near-IR spectrum done as in Section 4 yields that its blackbody temperature is T ∼ 1200 K (Fig. 4, bottom panel) and so slightly lower than the average value of the total dust emission of T ∼ 1300 K.
As an alternative explanation of the observed rms near-IR spectrum in Fig. 4, it is conceivable that the deficit in fractional variability at wavelengths ∼ > 1.2 µm is not due to the presence of a secondary, nonvariable dust component, but rather introduced by anticorrelation of the variability of the near-IR accretion disk flux and the hot dust flux, with the latter vary- The mean (black) and variable near-IR spectrum (red) for our campaign normalised at rest-frame 0.9 µm.The spectrum of a standard accretion disk (black dotted line) is traced to ∼ 1.2 µm in the variable component and only to ∼ 1 µm in the mean spectrum.Bottom panel: The spectral decomposition of the variable near-IR spectrum (black) into an accretion disk spectrum (black dotted line) and a hot dust component (red spectrum) yields a blackbody temperature of T ∼ 1200 K for the latter (red dotted curve).
ing in delayed response to the former.In such a case, we expect the total variability signal in the near-IR to be a superposition of the two variability signals and so the spectral shape of the rms spectrum to depend on the interplay between them.This superposition can then either enhance or weaken the total variance relative to the mean.In the following, we have simulated this situation with the DEMC formalism assuming that the secondary, variable near-IR component is the accretion disk itself.
The DEMC algorithm was developed by Schnülle et al. (2015) to suit multi-band optical/near-IR photometric dust reverberation campaigns and it can effectively construct the rms spectrum over this entire wavelength range.The analysis is carried out in two steps.First, the structure function parameters for the accretion disk variability are derived from modelling the optical lightcurve by maximum likelihood based on the method of interpolation and reconstruction of noisy, irregularly sampled data by Rybicki & Press (1992).The covariance function of the underlying Gaussian Process (GP) model is assumed to be a power-law (e.g.Press et al. 1992;Schmidt et al. 2010;Hernitschek et al. 2015) and the corresponding structure function can be written as: where A is the amplitude of the variability on a oneyear timescale and γ is the gradient of this variability.Secondly, a model consisting of the sum of a power-law, which represents the accretion disk emission, and a single blackbody, which represents the hot dust emission, is fit to the simultaneous multi-band optical/near-IR data.This model is: (3) where f g (t) and α are the variable accretion disk flux in the optical and the spectral power-law index, respectively, h is the Planck constant and k B is the Boltzmann constant.The dust model temperature T (t) evolves with time from an initial value T 0 as dT (t)/T (t) = ν • 0.25 df g (t − τ )/f g (t − τ ) in response to the variable accretion disk flux f g (t − τ ).The accretion disk flux is time-shifted by τ and processed with an efficiency ν (with 0 < ν < 1).The posterior distribution of the vector consisting of six parameters x = (C 1 , C 2 , ν, α, T 0 , τ ) is sampled by the Differential Evolution Markov Chain (DEMC; Ter Braak 2006) algorithm, which is basically a MCMC algorithm with multiple chains run in parallel.The algorithm evaluates the posterior probability density function for each iteration step and for each chain given uniform priors within sensible limits.
For the simulations, we assumed a fixed campaign length of 1000 days and modelled the input light-curve with stochastic short-term variability overlaid with a roughly sinusoidal long-term variability with a period length of ∼ 400 days.As the light-curves presented in Miller et al. (2022) show, this is a reasonable assumption for the variability pattern of the accretion disk in Mrk 876, which exhibits a quasi-oscillatory trend over a time period of about three years.Their data cover two peaks and two troughs of the long-term variability pattern, indicating roughly a period length of about twice the dust response time.Then, we chose response times of the hot dust emission relative to the variable (optical) accretion disk flux of τ = 50, 100, 150 and 200 light-days.In addition, we varied the processing efficiency of the dust and assumed values of ν = 0.5, 0.8 and 1, with the latter producing the maximum possible dust flux variability.Fig. 6 shows our results.It is clear that this model of two variable near-IR components can qualitatively reproduce the observations.The strongest depression of the hot dust variations is achieved for a combination of response time that gives maximal weakening of the total near-IR variability signal, i.e. half the period length and so ∼ 200 light-days, and an intrinsically reduced fractional variability of the hot dust emission due to a low processing efficiency.As the accretion disk flux contribution to the total near-IR flux decreases with wavelength (see Fig. 3), so does its weakening effect on the dust variability signal and the dust blackbody can always dominate the rms spectrum at the longest wavelengths.Then, if a secondary variable near-IR component is indeed present at ∼ 2 µm and it is the accretion disk, the implied outer radius is r out ∼ 4300 r g ∼ 56 light-days.Such an accretion disk is likely to harbour dust, which can contribute to the total blackbody emission.We will argue this case further in Section 6.3.Finally, we note that a reduced value of ν in our simulations has a similar effect as adding a constant dust component, i.e. it decreases the dust fractional variability.Therefore, it is conceivable that the observed rms spectrum of Mrk 876 is a combination of both scenarios described above.

THE DUST STRUCTURE IN MRK 876
It is of high interest to understand the chemical composition and grain size distribution of the circumnuclear dust in AGN, which can ultimately reveal how cosmic dust forms and evolves in different environments.If we can also constrain where the dust forms, we can assess its relationship to and influence on the other central structures of the AGN, such as, e.g., the BLR and accretion disk.As Landt et al. (2019) showed, a spectroscopic near-IR reverberation mapping campaign that uses cross-dispersed data covering a very wide wavelength range is ideally suited to investigate these questions.In the following, we discuss our results for Mrk 876 with a focus on how the combination of a luminosity-based dust radius with a contemporaneous dust response time can constrain the astrochemistry of the dust (Section 6.1), how the simultaneous measurement of the dust temperature and dust radius can shed light on the dust geometrical structure (Section 6.2), and the potential of a near-IR variable (rms) spectrum to unravel the structure of accretion disks in AGN (Section 6.3).Finally, we relate our results to AGN studies at other wavelengths and develop a new paradigm that connects AGN and protoplanetary discs (Section 6.4).

The astrochemistry of the hot dust
From the radiative equilibrium relationship (Eq. 1) it is clear that if estimates of the dust temperature T and irradiating luminosity L uv are available and the dust radius can be measured by an independent method, the dust emissivity parameter Q em , which strongly depends on the dust chemical species and grain size, can be constrained.Independent measures of the dust radius can come from, e.g., the dust response time measured through reverberation mapping, a geometric measurement based on near-IR interferometry or a polarimetric estimate of the location of the scattering region.
Our near-IR spectroscopic data does not extend over a time period long enough to be able to derive a reliable dust response time for Mrk 876.However, the WISE monitoring data for the years 2010-2018 analysed by Lyu et al. (2019) overlap with our campaign and the W 1 band (3.4 µm) is expected to sample the same hot dust component by covering its flux close to the SED maximum (see Fig. 3).Lyu et al. (2019) estimated a dust response time of τ = 297 ± 17 light-days (observer frame) for the W 1 band, assuming no corrections (see their Fig. 5).Adjusting the initial guesses for their fitting routine based on cross-correlation function results, the revised value was τ = 372 ± 63 light-days (see their Table 4).Their result is consistent with the distance to the equatorial scattering region of 254 ± 39 light-days (rest-frame) estimated by Afanasiev et al. (2019)  This independent measure of the dust radius of ∼ 300 light-days is smaller than all average dust radii we estimated from thermal equilibrium considerations.The least discrepancy is found relative to the case of a pure blackbody ( R d,lum = 469 ± 30 light-days), which corresponds to dust composed of large grains.The average luminosity-based radius for dust of small grains is much larger and exceeds the dust response time by a factor ∼ > 10 and ∼ > 16 for carbon and silicate dust, respectively.Therefore, it seems that the hot dust in Mrk 876 is dominated by grains with sizes of at least a few µm, as found for NGC 5548 by Landt et al. (2019).However, contrary to the case of NGC 5548, the luminosity irradiating the dust in Mrk 876 appears to be only a fraction of the bolometric luminosity that we see.Alternatively, if the dust is exposed to the total bolometric luminosity, something must reduce its temperature.In this case, it is worth noting that the response-weighted dust radius predicts a dust temperature of T ∼ 1700 K, which is close to the sublimation temperature for carbonaceous dust (T sub ∼ 1800 − 2000 K; Salpeter 1977).Carbonaceous dust was found to be the dominant species in the hot dust of NGC 5548.
It was previously postulated that the circumnuclear dust in AGN is dominated by large grains (e.g.Laor & Draine 1993;Maiolino et al. 2001).The spectroscopic near-IR monitoring campaigns on NGC 5548 and Mrk 876 have shown that this is indeed the case.Furthermore, they put an important lower limit on the dust grain size of a few µm, which suggests that the hot dust is assembled in a disk-like structure.When the dust grains are very small, their coupling to the ambient gas is strong and Brownian motion dominates the collision rate.Brownian-motion driven agglomeration can efficiently grow µm-sized dust particles, but the induced collision velocity decreases with increasing aggregate mass.Then, as dust grains grow they increasingly decouple from the gas and differential settling becomes the main driver for grain collisions and so further grain growth (Weidenschilling 1977;Blum & Wurm 2008).Therefore, large dust grains are expected to be found mainly in the midplane of a disk-like structure.If the gas is turbulent, we expect it to drive diffusion of dust in the vertical direction, which will oppose the effect of vertical settling and give the disk-like dust structure a scale height.

The enlarged dust-free inner hole
The lower dust reponse time in Mrk 876 relative to the dust radius estimated from radiative equilibrium considerations contrasts with the findings of Landt et al. (2019) for NGC 5548, where the two dust radii were similar.However, it is in line with general comparisons between interferometric dust radii and those determined by photometric reverberation campaigns; the former are found to be larger than the latter by a factor of ∼ 2 (Kishimoto et al. 2009(Kishimoto et al. , 2011;;Koshida et al. 2014).Since interferometry effectively measures the luminosity-weighted radius, this result can be understood if one assumes that response-weighted radii are dominated by dust located at the smallest radii, whereas luminosity-weighted radii are dominated by dust grains with the highest emissivity, i.e. the hottest and largest grains.Due to their lower heat capacity, the highest temperature large grains can attain will always be lower than the highest temperature small grains can be heated to, and so their emission will come from further out.We can easily rule out this scenario for Mrk 876 since the blackbody temperatures in the mean and rms near-IR spectrum are similar.
Alternatively, Kawaguchi & Mori (2010, 2011) ascribed the observed difference between response-and luminosity-weighted dust radii to a concave (or bowlshaped) dust geometry.Due to the anisotropy of the accretion disk emission, dust located at larger viewing angles relative to the irradiating UV/optical flux will see it reduced by a factor cos θ, with θ the angle between the accretion disk rotation axis and the location of the dust.The interferometric radii would then be dominated by the region containing the bulk of the dust mass, which in this geometry is located further out, whereas reverberation would simply sample the smallest and so most variable radii.In this scenario, the concave shape of the dust structure is caused by the assumption of a constant maximum dust temperature, taken to be the dust sublimation temperature, which is found further away the higher the accretion disk luminosity appears to the dust.Our observations for Mrk 876 can be well explained with this scenario.The variable and mean near-IR spectra yield similar dust temperatures, albeit below the sublimation temperature.Then, in order to make the luminosity-based dust radius consistent with the dust response time, we require that the variable dust sees the bolometric luminosity reduced by a factor of ∼ 4. If our viewing angle to the accretion disk is θ = 0 • , as assumed in our calculations of the accretion disk spectrum (Section 4), a variable dust component in a structure inclined at θ ∼ 75 • and so closely aligned with the plane of the accretion disk would fulfill this requirement.The large dust deficit observed in the rms spectrum relative to the mean spectrum (Section 5) would then be due to the fact that the total emission is dominated by the much larger dust mass located further out, which appears constant since it varies on much longer timescales than our campaign is sensitive to.
An indication that the accretion disk illumination anisotropy considered by Kawaguchi & Mori (2010, 2011) affects low-and high-luminosity AGN differently was also found by Minezaki et al. (2019).Their sample, which included NGC 5548 and Mrk 876 at the low-and high-luminosity ends, respectively, and spanned about four orders of magnitude in luminosity, showed a bestfit slope for the logarithmic reverberation dust radius versus optical luminosity relationship of ∼ 0.4, i.e. shallower than the slope of 0.5 predicted by thermal equilibrium considerations.Interestingly, using the dust response time and average accretion disk luminosity for NGC 5548 determined by Landt et al. (2019) and the values of these parameters for Mrk 876 from our study we derive also a slope of ∼ 0.4.Clearly, it would be worthwhile to further populate the undersampled top end of the relationship between dust radius and optical continuum luminosity in order to investigate how the dust irradiation pattern depends on accretion disk luminosity.
But one of our most important results cannot be explained solely with the scenario of Kawaguchi & Mori (2010, 2011), although the implied flared disk-like structure for the hot dust helps.Similarly to what Landt et al. (2019) observed in NGC 5548, the hot dust in Mrk 876 is not close to its sublimation temperature.This finding implies the presence of an enlarged dustfree inner region in AGN and so a luminosity-invariant inner edge of the torus, i.e. a "torus wall".Since the sublimation temperature for carbon dust is higher than for silicate dust (∼ 1900 K and ∼ 1400 K, respectively), the extent of this 'inner hole' would be considerably larger if indeed the former chemical species dominates the dust composition.Such a dusty wall was confirmed for NGC 4151 with several epochs of both interferometric and reverberation dust radius measurements (Koshida et al. 2009;Pott et al. 2010;Schnülle et al. 2013).This paradigm for the AGN dust structure is similar to the standard picture for protoplanetary disks around young stars.Interferometry finds in many of these disks that the inner dust radius is significantly larger than the dust sublimation radius.This cavity, which is commonly referred to as the 'inner hole', is believed to be filled with gaseous disk material that can change its optical thickness and thus influence the location of the puffed-up 'wall' or 'inner rim' of the dusty, flared and passively illuminated protoplanetary disk (Monnier et al. 2005;Dullemond & Monnier 2010).Our result that circumnuclear AGN dust is composed of relatively large grains was found also for protoplanetary disks (van Boekel et al. 2004;Kóspál et al. 2020).We will explore in more detail the similarity between these two classes of astrophysical objects in Section 6.4.

A secondary hot dust component
In Section 5, we showed that the observation of a reduced dust emission in the rms near-IR spectrum of Mrk 876 relative to its mean spectrum can also be explained by the accretion disk spectrum extending well into the near-IR.The implied large disk radii were ∼ 25 and ∼ 56 light-days in the J and K bands, respectively.These values are similar to estimates from a contemporaneous accretion disk reverberation mapping campaign.Miller et al. (2022) recently presented the accretion disk analysis using optical photometric monitoring in several bands (u, g, r, i, and z) conducted between 2016 Mar and 2019 May, including our data set presented in Section 3.2.They found that the response times depend on the wavelength as expected for a standard thin disk.Then, by extrapolating their optical response times to the near-IR, they estimated outer disk radii of ∼ 20 and ∼ 45 light-days in the J and K bands, respectively.Given these extents, it is of interest to understand if they are large enough to potentially harbour dust.
We can estimate the location of the accretion disk dust, if it exists, as follows.First, we used the simultaneous measurement of the dust temperature in the rms spectrum and the dust response time for the variable hot dust component together with the emissivity law we determined in Section 6.1 to "calibrate" the temperature-radius relationship for the dust and thus to estimate the bolometric luminosity that produces it.We show this relationship, which is that of a blackbody and so has the form T ∝ R −0.5 , in Fig. 7 (red solid line).As already discussed in Section 6.2, the corresponding total bolometric luminosity that seems to irradiate the dust is log L bol ∼ 45.5 erg s −1 , which is a factor of ∼ 4 lower than the average bolometric luminosity of log L uv = 46.09± 0.03 erg s −1 estimated from fits to the near-IR spectra.Assuming that this irradiating luminosity is produced by an accretion disk, we then calculated the temperature-radius relationship for it, which is of the standard form T ∝ R −0.75 and so steeper than that for the hot dust (black solid line).These curves immediately show that dust is expected to condense out in the accretion disk gas at much smaller radii than those possible for gas exposed to the full bolometric luminosity.Then, if this dust was driven off the accretion disk, e.g., by a dusty outflow, it would be quickly destroyed at considerable heights (Czerny et al. 2017).For the observed dust temperature, we expect the accretion disk of Mrk 876 to become dusty at radii ∼ > 40 light-days.This secondary hot dust will not necessarily be externally illuminated and so reverberate in the classical sense.Instead, if it is composed of µm-sized grains, it will be well-coupled to the gas and vary on the same timescales.Dust with a variability on such short time-scales is likely to appear constant in our rms near-IR spectrum.
An additional argument for the reality of a secondary hot dust component in Mrk 876 can also be derived from the total observed dust luminosity.The values listed in Table 2 (column (4)) give an average of log L d = 45.55 ± 0.02 erg s −1 , which, if interpreted as being produced by a single dust component (e.g. with the geometry proposed by Kawaguchi & Mori (2010)), implies a dust covering factor of unity given the "calibrated" bolometric luminosity irradiating it.Given that the hot dust is optically thick to the UV/optical accretion disk luminosity, such an AGN would be an obscured, Compton-thick one.On the other hand, if we assume an additional dust component in the mean spectrum, the rms spectrum tells us that the variable dust component makes up only roughly a third of the total observed dust flux and so has a covering factor of ∼ 30%.We note that by coincidence we would estimate a similar value for the covering factor also from the mean spectrum alone, since both the bolometric luminosity illuminating the torus dust and the dust luminosity would be overestimated by a similar factor.Then, if the dust deficit in the rms spectrum appears exacerbated due to the influence of the variable accretion disk (and possibly its associated dust), as discussed in Section 5, the torus dust covering factor would be higher.
The first tentative evidence for dust in the accretion disk was found for NGC 5548 by Landt et al. (2019).If it can be shown that such dust is ubiquituous in AGN, it would provide a possible connection to the BLR within models for radiatively accelerated dusty outflows launched from the outer accretion disk (e.g.Elitzur & Shlosman 2006;Czerny et al. 2016Czerny et al. , 2017;;Esser et al. 2019).In this case, carbon dust would provide the largest opacity and lead to an inflated disk structure (Baskin & Laor 2018).Support for this scenario in Mrk 876 comes from the remarkable agreement between the location we estimated for the accretion disk dust and the Hβ BLR radius of ∼ 50 light-days measured by Bao et al. (2022) for their reverberation mapping campaign conducted between 2016 Dec and 2018 May.Disk-like dust structures in AGN have started to be revealed also by interferometric observations achieving an unprecedented high spatial resolution in the mid-IR (Isbell et al. 2022).

Are AGN and protoplanetary disks similar?
Large accretion disks that extend well beyond the selfgravity radius are now routinely inferred for AGN by accretion disk reverberation studies (see review by Cackett et al. 2021).Evidence for the presence of dust in these disks would imply even greater radii, thus exacerbating the stability issue.A possible solution and an alternative explanation for our finding of a significantly reduced bolometric luminosity of the accretion disk in Mrk 876 as seen by the dust is that it is a passive disk.A passively illuminated disk, such as a protoplanetary disk, has a temperature-radius relationship indistinguishable from that of a standard accretion disk that produces its own radiation via gravitational energy release, i.e. it is of the same form T ∝ R −0.75 (e.g.Cackett et al. 2007).Accretion disk flux variability in AGN is readily interpreted as reprocessed emission from a central, highly variable, hard X-ray source since the observed time-scales are much shorter than the viscuous timescale, i.e. the time it takes for a change in accretion rate to trigger a change in disk temperature and so flux.The passively irradiated flux is usually assumed to be only a small fraction of ∼ 10 − 20% of the total disk emission (Starkey et al. 2016;Gardner & Done 2017).The value of ∼ 25% that we estimated for Mrk 876 is close to this range.An advantage of the passive disk scenario over the cos θ-effect in a bowl-like geometry, such as the one proposed by Kawaguchi & Mori (2010), is that it predicts the onset of dust in the accretion disk around the Hβ BLR radius; if only an anisotropic accretion disk emission as seen by the dust is assumed, the accretion disk itself would have the total (unreduced) bolometric luminosity and so be hotter, leading to an onset of dust much further out.We also note that the relative scaling we obtained between the mean and rms near-IR spectra in the optical regime of a factor of ∼ 6 (Section 5) is similar to the factor of ∼ 4 reduction in total bolometric luminosity required by the variable dust.Therefore, it could well be that the dust is irradiated only by the variable, reprocessed disk flux.
Passive disks can be distinguished from standard accretion disks, since, as a consequence of vertical hydrostatic equilibrium, they are flared, with the disk relatively thicker at larger radii.The vertically isothermal flared disk diverges from the flat disk solution at intermediate radii and approaches T ∝ R −3/7 at large radii, i.e. it is very similar to a single blackbody (e.g.Kenyon & Hartmann 1987).The dust is crucial for the disk thermodynamics (and for planet formation), but represents only ∼ 1% of the total mass and is easier to observe than the gas.Self-gravity is generally only significant in very massive disks and protoplanetary disks are often not massive enough to become gravitationally unstable.It is then conceivable that in Mrk 876 the sum of both temperature-radius profiles displayed in Fig. 7 (black and red solid lines) represents a single entity, namely, a passively illuminated, flared and dusty disk.A wind off this disk could then be the origin of the BLR.

SUMMARY AND CONCLUSIONS
We have conducted the first spectroscopic near-IR monitoring campaign on Mrk 876, one of the intrinsically most luminous AGN in the nearby Universe.Our cross-dispersed spectroscopy can measure dust temperatures with high precision, which allows us to derive the luminosity-based dust radius for different grain properties.When comparing it to an independent measure of the dust radius, we can then constrain the astrochemistry of the hot dust.Furthermore, we can construct the variable (rms) near-IR spectrum over a relatively large wavelength range.
Our main results can be summarised as follows.
(i) Assuming thermal equilibrium for optically thin dust, we find that the luminosity-based dust radii are larger than the dust response time obtained by a contemporaneous photometric reverberation mapping campaign, with the least discrepancy (of a factor of ∼ 2) found relative to the result for a wavelength-independent dust emissivity law, i.e. a blackbody, which is appropriate for grains of relatively large sizes (of a few µm).This result can be well explained by a flared, disk-like structure for the hot dust, whereby the anisotropy of the accretion disk emission causes a decrease in illumination (by a factor of ∼ 4 in our case), as first proposed by Kawaguchi & Mori (2010, 2011).
(ii) The near-IR variable (rms) spectrum tracks the accretion disk spectrum out to longer wavelengths (of ∼ 1.2 µm) than the mean spectrum (of ∼ 1 µm) due to a reduced dust emission in the former.The implied outer accretion disk radius is much larger than the self-gravity radius and consistent with the extrapolated results from the contemporaneous, multi-band optical accretion disk reverberation mapping campaign of Miller et al. (2022).The flux variability of the hot dust with respect to the accretion disk is reduced by a factor of ∼ 3 and could be due to either the presence of a secondary hot dust component in the mean spectrum or the destructive superposition of the dust and accretion disk variability signals or some combination of both.
(iii) The large extent of the accretion disk, which is likely to harbour dust in its outer regions, and the low bolometric luminosity as seen by the hot dust can also be explained if we assume that AGN disks are similar to protoplanetary disks around young stars.A passively illuminated, flared and dusty disk would naturally provide a single, continuous structure for the temperatureradius profile determined by AGN accretion disk reverberation studies (T ∝ R −3/4 ) at small radii and that of the hot dust blackbody (T ∝ R −1/2 ) determined here at large radii.
HL is indebted to Brad Peterson for suggesting the science target and leading the infrared observing proposals and to Peter Abraham, Agnes Kospal and Hubert Klahr for generously sharing their knowledge about protoplanetary disks.HL and JAJM thank Michael Fausnaugh for help with applying his mapspec routine to the near-IR spectroscopy.HL acknowledges a Daphne Jackson Fellowship sponsored by the Science and Technology Facilities Council (STFC), UK.J.A.J.M acknowledges the support of STFC studentship ST/S50536/1.HL, JAJM and MJW acknowledge support from STFC grants ST/P000541/1 and ST/T000244/1.KH and JVHS acknowledge support from STFC grant ST/R000824/1.ERC acknowledges the support of the South African National Research Foundation.
Facilities: Gemini(GNIRS), LCOGT(optical) Figure 1.Gemini GNIRS near-IR spectrum from 2016 May 25, shown as observed flux versus rest-frame wavelength.Emission lines listed in Table4ofLandt et al. (2008) are marked by dotted lines and labeled; black: permitted transitions, green: permitted Fe II multiplets (not labeled), red: forbidden transitions and cyan: forbidden transitions of iron (those of [Fe II] not labeled).

Figure 2 .
Figure 2. Top two panels: LCOGT g and zs band light-curves during the period 2016/17.Bottom two panels: Gemini GNIRS near-IR spectral light-curve around the observed wavelength of 8700 Å (black filled circles), both original and corrected using PrepSpec photometric correction factors based on the [S III] λ9531 line, versus the zs band light-curve (red circles).The zs band light-curve corrected for a constant host galaxy contribution to match the near-IR spectral fluxes is also shown (green filled circles).
Figure 3.Gemini GNIRS near-IR spectrum from 2016 May 25 shown as luminosity versus rest-frame wavelength (black spectrum).We decomposed the continuum into two components, namely, an accretion disk spectrum that approximates the wavelength range ∼ < 1 µm (black dotted curve) and still dominates at 1 µm (vertical dashed line) and hot dust emission (red spectrum).We fitted the hot dust continuum with a blackbody spectrum (red dotted curve) and modified blackbody spectra for carbon and silicate dust (green and blue dotted curves, respectively) with resulting best-fit temperatures as listed in Table2.

Figure 4 .
Figure 4. Top panel:The mean (black) and variable near-IR spectrum (red) for our campaign normalised at rest-frame 0.9 µm.The spectrum of a standard accretion disk (black dotted line) is traced to ∼ 1.2 µm in the variable component and only to ∼ 1 µm in the mean spectrum.Bottom panel: The spectral decomposition of the variable near-IR spectrum (black) into an accretion disk spectrum (black dotted line) and a hot dust component (red spectrum) yields a blackbody temperature of T ∼ 1200 K for the latter (red dotted curve).

Figure 5 .
Figure5.Top panel: The flux difference between the normalised mean and variable spectra (blue) can be attributed to an additional (non-variable) dust component with a blackbody temperature of T ∼ 1300 K (blue dotted line).The template of an elliptical host galaxy is also shown (magenta).Bottom panel: The fractional amplitude of the variations as a percentage of the mean flux versus the rest-frame wavelength.

Figure 6 .
Figure 6.Simulations of the mean (black line) and rms spectra (coloured lines) with the DEMC formalism assuming a secondary variable component at near-IR wavelengths, which is produced by the accretion disk emission.The rms spectra are simulated for four different dust response times (τ = 50, 100, 150 and 200 light-days, indicated by the green, red, light blue and dark blue solid lines, respectively) and three different values for the dust processing efficiency (ν = 1, 0.8 and 0.5, shown in the left, middle and right panel, respectively).
based on spectropolarimetric observations from 2015.Furthermore, the dust response time in Mrk 876 appears to remain relatively constant.From an optical/near-IR photometry campaign during the years 2003-2007, Minezaki et al. (2019) obtained values similar to Lyu et al. (2019) of τ = 320 ± 34, 327 +42 −36 and 327 +41 −35 light-days (observer frame), assuming three different powerlaw spectral slopes for the accretion disk spectrum when decomposing the flux in the K band.

Figure 7 .
Figure 7.The temperature-radius relationships for the hot dust (red curve) and the accretion disk illuminating it (black curve).They correspond to a bolometric luminosity of log L bol = 45.5 erg s −1 , which we derived from "calibrating" the relationship for the hot dust to its observed temperature in the rms near-IR spectrum and the rest-frame dust response time (T = 1200 K and R d,rev = 300 light-days, shown as the red horizontal and vertical dashed lines, respectively).The blue dashed lines indicate the accretion disk reverberation radii of ∼ 20 and ∼ 45 light-days estimated by Miller et al. (2022) for the J and K bands, respectively.The black vertical dashed line indicates the Hβ BLR response time of ∼ 50 light-days measured by Bao et al. (2022).

Table 1 .
Gemini journal of observations

Table 2 .
Physical parameters for the calculation of luminosity-based dust radii