The turbulent variability of accretion discs observed at high energies

We use numerical stochastic-viscous hydrodynamic simulations and new analytical results from thin disc theory to probe the turbulent variability of accretion flows, as observed at high energies. We show that the act of observing accretion discs in the Wien tail exponentially enhances small-scale temperature variability in the flow, which in a real disc will be driven by magnetohydrodynamic turbulence, to large amplitude luminosity fluctuations (as predicted analytically). In particular, we demonstrate that discs with more spatially coherent turbulence (as might be expected of thicker discs), and relativistic discs observed at larger inclinations, show significantly enhancement in their Wien-tail variability. We believe this is the first analysis of relativistic viewing-angle effects on turbulent variability in the literature. Using these results we argue that tidal disruption events represent particularly interesting systems with which to study accretion flow variability, and may in fact be the best astrophysical probes of small scale disc turbulence. This is a result of a typical tidal disruption event disc being naturally observed in the Wien-tail and likely having a somewhat thicker disc and cleaner X-ray spectrum than other sources. We argue for dedicated X-ray observational campaigns of tidal disruption events, with the aim of studying accretion flow variability.


INTRODUCTION
Accretion, whereby the energy of infalling material is liberated by turbulent dissipation and ultimately radiated away, is one of the most ubiquitous and important phenomena in the universe.The central object powering the accretion process can range from young stellar objects (YSOs), white dwarfs (WDs) to black holes (BHs).Similarly, there are many different possible mechanisms for the origin of the accreted material, including Roche-lobe overflow in cataclysmic variables (CVs) and X-ray binaries (XRBs), tidal disruption events (TDEs) where a star on a close orbit of a BH is tidally disrupted, and the accretion of galactic gas by supermassive BHs (SMBHs) that power active galactic nuclei (AGN).Understanding the nature of the accretion process is important not only for understanding these objects, but also because of the role accretion has in wider astrophysics, for example in the triggering of Type Ia supernovae and in the effect feedback from AGN has on the evolution of their host galaxies.
One interesting property of these systems is that, despite the large range of scales that they cover, they show remarkably similar variability properties.These include a log-normal luminosity distribution; a linear relation between the mean flux varying over long timescales and the root mean square (rms) variability on short timescales; broad noise across a range of frequencies as seen in the power spectral density (PSD); and hard and soft lags observed between different energy bands.The presence of these features has been observed in YSOs (Scaringi et al. 2015), CVs (Scaringi et al. 2012;Veresvarska & Scaringi 2023), XRBs (Nowak 2000;Gleissner et al. 2004;Gandhi 2009) and AGN (Lyutyi & Oknyanskii 1987; ⋆ E-mail: andrew.mummery@physics.ox.ac.ukGaskell 2004;Arévalo et al. 2006;Markowitz et al. 2007;Fabian et al. 2009;de Marco et al. 2011;Vaughan et al. 2011).
The theory of propagating fluctuations, first proposed by Lyubarskii (1997), is often used to explain these observed phenomena.In this theory, local fluctuations in the accretion rate are created at all radii.These fluctuations then propagate inwards 1 , combining with those from smaller radii.Assuming that the timescale these fluctuations are created on increases at larger radii in the disc (note that, in a Newtonian disc, all the key timescales scale as r 3/2 where r is the radial distance from the central object), this naturally creates a broad PSD across a range of frequencies.Further assuming that the fluctuations combine multiplicatively, this also gives rise to a linear rms-flux relation, as the fractional variation created by the high frequency fluctuations in the inner disc provide a constant fractional modulation to the low frequency fluctuations from the outer disc.This linear rms-flux relation, provided that it holds for all frequencies, implies that the underlying lightcurve must be log-normally distributed (Uttley et al. 2005).As fluctuations propagate inwards, they will first pass through the cooler, outer regions of the disc, before the hotter, inner regions.This will naturally create hard lags, as any signal will appear first in the lower energy bands.Mushtukov et al. (2018) additionally showed that, if the theory is extended to allow for any fluctuations to propagate outwards, then this can also lead to soft lags.All of these results also hold in a relativistic theory of propagating fluctuations (Mummery 2023).
Since the seminal work of Lyubarskii (1997), there has been extensive analytic work on propagating fluctuations (e.g.Kotov et al. 2001;Arévalo & Uttley 2006;Ingram & Done 2011;Ingram & van der Klis 2013;Mushtukov et al. 2018;Mummery 2023).Additionally, numerical studies have allowed the model to be tested in the non-linear regime.These have been performed in 1D (based on the standard disc diffusion equation with a stochastic α parameter) (Cowperthwaite & Reynolds 2014;Turner & Reynolds 2021), and have also recently been generalised to 2D (Turner & Reynolds 2023).Additionally, high resolution magnetohydronamic (MHD) simulations have been found to exhibit behaviour which is consistent with the propagating fluctuations picture (Hogg & Reynolds 2016;Bollimpalli et al. 2020).
One of the more recent theoretical disc variability predictions (not precisely reliant on the theory of propagating fluctuations, merely the log-normality of disc temperature variability) was derived by Mummery & Balbus (2022), who highlighted a key result from thin disc theory.Mummery & Balbus (2022) argued that when observed in the Wien tail (i.e., at photon energies larger than the peak temperature of the disc), fluctuations in the disc temperature are exponentially amplified into significantly larger luminosity fluctuations.This is formally relevant for all disc systems, but in particular for discs formed in tidal disruption events, which typically peak at temperature kTp ≲ 100 eV, a factor of a few below typical X-ray bands (with lower bandpass energies typically E l = 300 eV).As these systems are observed in their Wien-tail, and therefore may well have exponentially enhanced variability, these could be the ideal disc systems to observe and probe our theoretical descriptions of disc variability in an astronomical setting.
Before comparing to observations however, it is imperative to ascertain a more precise set of predictions of the modern theory of propagating fluctuations (in particular the numerical simulations of Turner & Reynolds 2023), when filtered through the Wien-tail modelling of Mummery & Balbus (2022).Deriving observational predictions is the purpose of this paper.We shall do this by taking the numerical models of Turner & Reynolds (2023), specialising to regions of parameter space relevant for observations of tidal disruption events.We then "observe" these systems, demonstrating that the analytic theory of Mummery & Balbus (2022) reproduces the numerical results well, allowing us to calibrate observational predictions.
We will demonstrate that, if the canonical theory of propagating fluctuations is correct, that tidal disruption events should be observed to show large, and potentially up to an order of magnitude, fluctuations of their X-ray luminosity on timescales as short as days-weeks (or even potentially hours for lower mass black holes).This means that tidal disruption events represent ideal probes of the small-scale disc turbulence in black hole accretion discs.
The layout of the paper is as follows.In section 2 we recap the work of Turner & Reynolds (2023), and describe the numerical twodimensional viscous-hydrodynamic model used throughout this paper.In section 3 we recap the analytical model of Mummery & Balbus (2022), and show that the analytical model accurately reproduces the full numerical simulations.We then demonstrate that different assumptions regarding the properties of the stochastic viscosity (in particular the length scale over which the viscosity is coherent) produce pronounced effects on the Wien-tail lightcurves, and in particular we argue that thicker discs will produce significantly larger variability.In section 4 we highlight relativistic viewing angle effects, and show that at fixed observed X-ray luminosity more inclined relativistic disc systems will be observed to be more variable than face on systems (at the factor 2 level).In section 5 we discuss the relevance of our results for observations of tidal disruption events, before concluding in section 6.

NUMERICAL MODEL
The numerical model analysed in this work is presented in detail in Turner & Reynolds (2023), where a two-dimensional (vertically integrated) generalisation of the one-dimensional model presented in Turner & Reynolds (2021) was derived and analysed.We only briefly summarise the fundamental equations solved by the model here, and direct the interested reader to these papers for further information and discussion.
The disc fluid is assumed to evolve under the equations of viscous hydrodynamics, with a Newtonian gravitational potential.Explicitly, the governing equations are mass conservation (written in Cartesian index notation, where repeated indices are summed over) and momentum conservation In these expressions ρ is the disc density, ⃗ v the disc velocity, p the disc pressure and δij is the Kronecker delta.We use the notation xi for the Cartesian coordinate system and t is time as usual.The viscous stress tensor is given by and we assume a standard central point source gravitational potential where G is the Gravitational constant, M the mass of the central BH and R the spherical radius.As is now standard, the coefficient of disc viscosity is set by the Shakura & Sunyaev (1973) α-prescription, namely where cs = p/ρ and H are the disc's speed of sound and pressure scale height respectively.The work of Turner & Reynolds (2023) further assumed that the disc is sufficiently thin (i.e.H ≪ r, where r is the cylindrical radius in the plane of the disc) that the exact details of the vertical structure can be ignored, allowing the problem to be confined to a 2D, vertically integrated approach.This allows the density ρ to be replaced with the surface density Σ. Formally, this is done (in cylindrical polar coordinates (r, ϕ, z) which each have their standard meaning) as but this can be thought of simply as Σ ≡ ρH.A similar equation defines the vertically integrated pressure, P ≡ pH.It is implicitly assumed that other variables, such as the fluid velocity, are not a function of z.The thin disc assumption also allows R to be replaced with r in eq. 4.
In order to close the equations, the (approximate) solution of vertical hydrostatic equilibrium is used to relate the sound speed to the disc scale height as The disc aspect ratio H ≡ H/r is taken as an input parameter of the model, which also sets the sound speed.Note that the sound speed varies as r −1/2 through the disc.The simulations analysed in this work have an aspect ratio of H = 0.1, although other aspect ratios were explored by Turner & Reynolds (2023).Unlike traditional approaches, α is not taken to be constant in this analysis.Instead, in an attempt to mimic the effects of the real magnetohydrodynamic turbulence in the disc, the α-parameter is taken to be a stochastic random variable, and is equal to where a value of α0 = 0.1 is used in all the simulations of Turner & Reynolds (2023).β is a (two-dimensional) stochastic random variable, which is evolved stochastically in time and is additionally advected with the flow, in an effort to mimic the properties of the magneto-rotational instability (Balbus & Hawley 1991), and the subsequent turbulent evolution of the disc.To be precise, β is evolved according to an Ornstein-Uhlenbeck (hereafter OU) process.An OU process can be thought of as a generalisation of a random walk to continuous time, in which in addition the properties of the process have been chosen so that there is a tendency for the parameter to move back towards a mean value, with a growing attraction to the mean as the distance between the parameter and its mean increases.At a given point (r, ϕ) in the disc, the instantaneous stochastic evolution under an OU process is given by where 1/ω is the characteristic timescale of the OU process, µ is the mean value of β, dW ∼ N (0, dt) is a random normal variable and ξ is a constant which determines the magnitude of the variation, which can be written in terms of the rms-amplitude of β as ξ = 2ω ⟨β 2 ⟩.
All the simulations presented in Turner & Reynolds (2023) took µ = 0 and β 2 = 1 2 .Note that the choice of µ = 0 does not mean that ⟨α⟩ = α0 since the exponential function in eq. 8 has a positive skew.Multiple driving timescales/frequencies were tested in Turner & Reynolds (2023) (and an even wider range tested in the 1D work of Turner & Reynolds (2021)), but the simulations used in this work use the fiducial driving frequency of This corresponds to a timescale which is 10 (1/α0) times longer than the orbital timescale, and was chosen to be similar to the MRI dynamo timescale which appears to be the most important timescale for the evolution of the effective α parameter in MHD simulations (e.g.Hogg & Reynolds 2016).
In addition to the required normality of the random variable dW , the exact form was chosen such that this random noise was spatially coherent over a length scale lc within the disc.Motivated by the nature of a mean-field theory "turbulent viscosity", we might expect that (to first order) turbulent eddies within the disc will have a size given approximately by the disc thickness.Therefore, the fiducial assumption used in Turner & Reynolds (2023) is that lc = H, although other simulations where lc is some multiple of H were performed and are used in this work.For the full details of the nature of the dW and the general stochastic prescription, we directed the interested reader to Appendix A within Turner & Reynolds (2023).
In addition to being a stochastic variable, the parameter β is assumed to be advected with the flow.In other words, the existing β field is evolved with time through Under eq.11, existing structures within the β field are evolved, and regions of high (or low) β will be sheared out due to the presence of differential rotation within the disc.Eqs. 9 and 11 completely describe the evolution of β, and ultimately the local viscous stress tensor Πij of the flow.
The coupled equations (1-11) are solved by the code PLUTO Mignone et al. (2007).Note that while these equations have been presented in Cartesian notation, they are solved in polar coordinates (r, ϕ).The initial condition of the disc surface density is that of the analytical one-dimensional steady state solution Pringle (1981) where r⋆ = 6GM/c 2 is the inner edge of the disc, set by the innermost stable circular orbit (ISCO) of a Schwarzschild BH.Additionally, the velocity field was initialised with v ϕ equal to the local Keplerian velocity v ϕ = GM/r with vr = 0 and β = 0 everywhere.The fluid is first evolved during a run-in phase in which the stochasticity is turned on (equations 9 and 11).This phase lasts for a total length of time t = 0.4 × 10 6 tg, where tg ≡ GM/c 3 .This ensures that the β field will have settled into a statistically steady state before we start taking results.The disc then continues to evolve forward for a further t = 1.6 × 10 6 tg, where results are taken.We save the full state of the inner disc, defined by 6 < rc 2 /GM < 1000, each ∆t = 1000tg.We therefore have 1600 snapshots of the disc for which to analyse.

Emergent disc spectrum
In cylindrical coordinates with no z dependence, the local energy liberated from viscous dissipation D(r, ϕ, t) is equal to where the factor of 1/2 comes from the two (upper and lower) surfaces of the disc.Assuming that this dissipated energy is radiated locally, this defines the usual effective temperature T (r, ϕ, t) of the disc by where σ is the usual Stefan-Boltzmann constant.An example snapshot of temperature of the inner disc is shown in Fig. 1.The temperature in this Figure is re-scaled to a system corresponding to a black hole of mass M = 2 × 10 6 M⊙, and average accretion rate Ṁ0 = 6 × 10 20 kg/s (PLUTO solves the equations in a dimensionless unit system, which can then always be rescaled by appropriate choice of M and Ṁ0; Turner & Reynolds 2023).
Clearly it can be seen in Fig. 1 that the disc temperature is a stochastic random variable, with small length scale coherence to the temperature field (in this simulation the coherence length is equal to the disc scale height lc = H = 0.1r).While the average (over azimuth) temperature profile approximately follows what would be expected from one-dimensional thin disc theory, the local structure of the disc is much more complicated.We shall go on to demonstrate that at low observing energies this two-dimensional structure Figure 1.An example snapshot of the two dimensional (x and y) disc temperature profile.The central mass was taken to be M = 2 × 10 6 M ⊙ , and the average mass accretion rate Ṁ0 = 6 × 10 20 kg/s.It is clear that the disc temperature is a two dimensional stochastic random field.There is small length scale coherence to the temperature field, which in this case is of order the disc scale height lc ∼ H = 0.1r.The spectrum and lightcurves which result from observing this stochastic temperature field inherit the temperature's underlying randomness.
effectively averages out, but that this is not true at high observing frequencies (compared to the peak disc temperature).
In a purely Newtonian theory of gravity (i.e., one without gravitational or Doppler energy shifting of photons, which are assumed to travel on straight lines), the emergent spectrum (energy received per unit time per unit area per unit photon energy) can be written where h and k are the Planck and Boltzmann constants respectively, D is the distance to the observer, and i is the inclination of the disc as observed by the distant observer.The spectral luminosity is defined by LE ≡ 4πD 2 FE, and is independent of the observers distance.
For each snapshot of the disc temperature T (r, ϕ, t), the above integral can be performed.In Fig. 2 we show an example series of observed spectra (plotted in terms of the spectral luminosity ELE), as a function of observing energy E. Each curve corresponds to one of the 1600 disc snapshots recorded in the simulations of Turner & Reynolds (2023).These spectral snapshots include the relativistic effects we discuss in section 4, and so represent our most sophisticated model.The physical parameters of the system are chosen so that the system corresponds to a typical tidal disruption event system M = 2 × 10 6 M⊙, Ṁ0 = 6 × 10 20 kg/s, with an average 0.3 − 10 keV X-ray luminosity of L ∼ few × 10 41 erg/s.Clearly, the stochastic nature of the disc temperature results in a stochastic disc spectrum, and this stochasticity is an inherent function of observing energy.Namely, at low observing energies the two-dimensional stochastic structure of the disc temperature effectively averages out over the 2π azimuthal range of the disc, and the disc temperature is reasonably well described by the canonical one-dimensional thin disc theory (contrast the numerical spectra to the black-dashed classical theory prediction).However, at high observing frequencies (relative to the peak disc temperature), the spectrum is much more sensitive to the local  1).The black hole mass was set to M = 2 × 10 6 M ⊙ , and the average mass accretion rate Ṁ0 = 6 × 10 20 kg/s, so as to produce typical 0.3 − 10 keV X-ray luminosities typical of a tidal disruption event (a particularly relevant astrophysical system).It is clear that the fractional variability in the disc luminosity increases as a function of observing energy, and is particularly pronounced in the Wien tail.The conventional 1D disc theory spectrum for these parameters is shown by a black dashed curve.Note that for the purposes of computing light curves the disc spectrum is evaluated above 1 keV, but this is not shown here.
physics of the inner accretion flow, and the fractional variability in the spectrum is much larger.This is especially true in the Wien tail of the spectrum, where order of magnitude spectral changes can be observed (as predicted analytically; Mummery & Balbus 2022).
In addition, the luminosity observed across a band, defined by can also be computed.In Fig. 3 we show the 0.3 − 10 keV X-ray lightcurve computed by integrating the spectra in Fig. 2 across the X-ray bandpass (i.e., E l = 0.3 keV, E h = 10 keV).The lightcurve is notable for its large intrinsic variability, showing order of magnitude fluctuations in X-ray luminosity on timescales as short as the time between saved snapshots ∆t = 10 3 tg.Note that for typical tidal disruption event parameters the gravitational timescale is short tg = GM/c 3 ≃ 5(M/10 6 M⊙) seconds.It is therefore likely that short (∼ hour-day) timescale X-ray variability will be a feature of these events (we will discuss current observational evidence for such behaviour in later sections of this paper).
It is clear that the variability of accretion flows at high photon energies has the potential to be leveraged as a probe of small length scale variability in accretion flows, particularly for tidal disruption event discs which fortuitously fall into this observational regime.In the following section we introduce the analytical model of Mummery & Balbus (2022), which seeks to describe the probability density function of the high photon energy disc lightcurves.The purpose being to describe the variance seen in these systems, and link the observables to the local structure of the disc through a simpler (to analyse and compute) analytical model based on classical disc theory.

The high-energy luminosity probability density function
As can be clearly seen in Fig. 2, the fractional variability observed in the spectrum of an accretion disc is an increasing function of observing energy.The large-variance properties of the Wien-tail was predicted analytically by Mummery & Balbus (2022), and can be intuitively understood as the exponential enhancement of disc temperature fluctuations (e.g., Fig. 1) from the L ∼ exp(−E/kT ) spectral shape of the Wien tail.This intuitive explanation can be placed on a firmer theoretical footing by treating the disc temperature T (r, ϕ, t) as a random variable, and exploring the implications of this variability for the spectrum (eq.15) and disc luminosity (eq.16).Both of these variables now must of course be considered random variables (a randomness they inherit from the underlying turbulent disc properties).The analysis proceeds by deriving analytical approximations for the two spectral integrals (eqs.15, 16), before folding through the randomness of the temperature field.
In Mummery & Balbus (2021) it was shown that equations ( 15) and ( 16) can be evaluated by performing a Laplace expansion about the hottest region in the disc.The physical insight here is that only the very hottest disc regions contribute to the spectrum (or luminosity) in the Wien-tail, as the flux from all other disc regions is exponentially suppressed (see Balbus 2014, for further discussion).
Performing the Laplace expansion and then integrating results in the following expression for the disc luminosity (eq.16) Here we have defined Tp the hottest temperature in the accretion disc.The radius Rp corresponds to the radial location of the peak temperature, while the constant η depends on the inclination angle of the disc and the disc's inner boundary condition, and is limited to the range 3/2 ≤ η ≤ 5/2.The constant χ1 ≃ 2.19 for classical thin disc theory, although this assumes an axi-symmetric one-dimensional disc temperature profile.As this is not the case in our numerical simulations, the parameter Rp may differ by an order unity factor from the true radial scale of the disc temperature.This result clearly highlights an important effect of observing systems in the Wien-tail: only a small region of the disc contributes to leading order (the hottest disc temperature is the only disc variable which appears).This means that local turbulence will play an oversized role in Wien-tail luminosity variance, as the observable is sensitive to the small scale structure of the disc.While this may make such observations difficult to fit using traditional approaches, it does open up the possibility of probing small scale disc turbulence with observations.
It should be noted that equation ( 17) was derived in Mummery & Balbus (2021) assuming a temperature profile which depend on only one (radial) spatial dimensional (e.g., T (r, t)).Such an approximation is not required, and identical temperature scaling relationships (i.e., L ∝ T η exp(−E/kT )) can be derived for two dimensional temperature profiles (as shown in Balbus 2014; Mummery & Balbus 2020), with differing (but largely physically irrelevant) prefactors in this two-dimensional case.
The mathematical question is then, if the high energy luminosity is described by a (dimensionless) function of the form and X is a random variable, how is Y distributed?In the above expression Y = L/L0 is a normalised luminosity, X = kTp/E l a normalised temperature, and L0 is a disc temperature independent constant defined by The parameter E l is the lower bandpass energy from equation ( 16).
The answer to this question is (Mummery & Balbus 2022) where pX is the probability density function from which the random variable X (in this case the disc temperature) is drawn and W is the Lambert W function (Corless et al. 1996), which is defined as the solution w of the following equation and z ≡ η −1 y −1/η .the probability density function of temperature fluctuations followed a log-normal profile (a result with much observational support, e.g., Uttley et al. 2005).The log-normal distribution has the following form where this function is valid only for x > 0 (positive temperatures), and µN and σ 2 N are the mean and variance of the underlying normal distribution (these will be related to the mean µT and variance σ 2 T of the temperature distribution shortly).Under this assumption the X-ray luminosity should be distributed according to This model may now be tested with the viscous hydrodynamic disc models of Turner & Reynolds (2023) or with observations of astrophysical systems.If temperature fluctuations in astrophysical systems are not log-normally distributed then eq. ( 24) will not accurately describe the luminosity variability, but the more general result (equation 20) will still hold.The modelling of different temperature distributions in equation ( 20) could be leveraged as a future test of one of the key predictions of the theory of propagating fluctuations (the log-normality of emission).
Recall that the true luminosity variable L is related to the variable Y by a simple scaling L = L0Y , and that L0 depends on the system parameters only through the black hole mass M (which sets the radial scale Rp in equation [19]).The high energy luminosity distribution is therefore a four-parameter distribution which may be found by combining equation ( 24) with the definition of L0 (equation [19]): pL(l; Rp, η, µT , σT ) = L −1 0 pY (y = l/L0; η, µN , σN ), (25) where the mean (µN ) and variance (σ 2 N ) of the underlying normal distribution which appear in equation ( 24), are related to the mean (µT ) and variance (σ 2 T ) of the temperature distribution by

Modelling the disc variability
The question we pose in this section is whether the emergent highenergy luminosity of the stochastic-viscosity hydrodynamic disc model of Turner & Reynolds (2023) (specifically the model with lc = H = 0.1r) can be described by the analytical results of Mummery & Balbus (2022).
The luminosity probability density function (eq.24) is a function of four free parameters Θ = (η, µN , σN , L0).For a given X-ray lightcurve (e.g., Fig. 3), the likelihood of observing that lightcurve, given the set of free parameters of the theory, is given by definition by the expression where {lj} is the list of N luminosity values.This explicitly assumes that each observed luminosity is independent of all previous 0 .97 0 .9 9 1 .0 Figure 5.The correlations between the parameters of the X-ray luminosity probability density function (eq.24), and their posterior distributions, for the fit to the numerical distributions shown in Fig. 4. The fractional temperature variance σ T /µ T is of order ∼ 1.5%, despite the X-ray lightcurves showing significant (factor ∼ few) variability.
luminosities; an assumption which is only true if the luminosities are observed at times separated by a period longer than the typical driving timescale of the turbulence.For the models considered here, the driving timescale in the (relevant) inner regions is of order ∼ 150tg (see eq. 10), whereas we "observe" on timescales of 1000tg, and this approximation is therefore valid.
To test the ability of the analytical result (eq.24) to reproduce the numerical light curves produced by the model of Turner & Reynolds (2023), we generate X-ray lightcurves L(t) for a disc model with M = 10 6 M⊙ and Ṁ0 = 6 × 10 21 kg/s, corresponding to typical tidal disruption event values (with an average 0.3 − 10 keV X-ray luminosity of order LX ∼ 10 42 erg/s).Note that the actual values of the parameters here are unimportant, owing to the dimensionless unit system that PLUTO employs.As the two dimensional temperature structure is the same for each run (with a given set of turbulence parameters) shifting the temperature by will result in identical variability structure (about a different mean amplitude) when observed at an energy We generate lightcurves for a range of values of E l (but fixed disc parameters), chosen so that the typical average luminosities of the different bands span a factor of 100.The smallest value of E l is 0.5 keV, and the larger E l values are all then equally spaced up to E l = 1 keV.
To model these lightcurves, we fit the temperature parameters µT and σT , and the disc size parameter Rp by maximising the loglikelihood In Fig. 4 we display the numerical (coloured histograms) distribution of disc luminosities (e.g., histograms of the disc lightcurves like those displayed in Fig. 3).As purple solid curves we display a fit to the analytical probability density function (eq.24), with one set of disc parameters (µT , σT , Rp) fit to all of the distributions simultaneously.We fix η = 5/2, which was found to best describe the numerical data.Physically this value of η corresponds to a temperature profile with sharp gradients, which is perhaps unsurprising given the temperature structure seen visually (Fig. 1), see Mummery & Balbus (2021) for further discussion on the physics of η.It is clear that the analytical distribution (eq.24) well describes the numerical results.This is an important result.It is not obvious that the variable luminosity resulting from a two-dimensional viscous-hydrodynamics simulation, with a stochastic viscosity, would be well described by simple analytic disc theory.As can be seen in eq. 13 the disc temperature is a non-linear function of viscosity, density and the velocity gradients of the flow, all of which are implicitly dependent on the stochastic variable β.The fact that this complex behaviour can be quite accurately described by thin disc theory coupled to a lognormal temperature profile is non-trivial, particularly as the local accretion rate and local temperature are no longer tightly correlated in these simulations (see Figure 20 of Turner & Reynolds 2023, for more details).The fact that this simple description works well opens up the possibility of using the theory (eq.24) to probe the properties of accretion disc turbulence, or alternatively help explain some of the variability of those sources observed in the relevant limit of E l ≫ kTp.
The values of (and correlations between) the best fitting parameters of the luminosity distribution are displayed in Fig. 5.While fitting we normalise the radial coordinate by rg = GM/c 2 , where M is the (exact) black hole mass variable used in generating the lightcurves, so that Rp represents a typical (dimensionless) length scale of the peak disc temperature.An extremely important result is Lognorm(µ T , σ T ) Figure 6.The distributions of the fixed azimuth maximum temperatures Tmax, and the angle averaged maximum temperature ⟨Tmax(t)⟩ ϕ (produced with M = 10 6 M ⊙ , Ṁ0 = 6 × 10 21 kg/s, see text).The small scale temperature variability (i.e., the one probed by Tmax(ϕ, t)) is significantly more variable than the angle-averaged temperature scale ⟨Tmax(t)⟩ ϕ .However, as can be seen by comparing the fitted distribution (black dashed curves) to the physical distributions, the light curve variability is (at least in a Newtonian system) dominated by the variability of the angle-averaged temperature profile.
that the fractional variability of the disc temperature is only of the percent level (σT /µT ≃ 1.5%), despite the X-ray variability of the light curve being at the factor ∼ 2 level for lower observing frequencies (blue histogram Fig. 4), and the factor ∼ 5 level for higher observing frequencies (red histogram, Fig. 4).This highlights the exponential amplification of disc temperature variability when observing in the Wien-tail.
It is interesting to compare the fitted distribution of the peak disc temperature (i.e., a log-normal distribution with mean µT and standard deviation σT ) with the physical properties of the disc system.Clearly (e.g., Fig. 1) the disc temperature is a complicated twodimensional random field.We can define two different "peak disc temperatures", the first being the maximum temperature in the disc at a given azimuth at a fixed time (the rationale here being that the angular structure probes the turbulent variation, as the radial temperature dependence is driven principally by the gravitational energy liberation) while we can also define an angle-averaged maximum temperature which is likely more closely related to the classical disc theory prediction.We plot the two distributions in Fig. 6, which highlights the difference in these two values and the structure of the accretion flow.The small scale temperature variability (i.e., the one probed by Tmax(ϕ, t)) is significantly more variable than the angle-averaged temperature scale ⟨Tmax(t)⟩ ϕ .Indeed, the small scale temperature structure varies by a factor of ∼ 2 over the entire simulation, while the averaged temperature shows significantly smaller variability.Also shown in Fig. 6 is the best fitting temperature distribution fit to the X-ray lightcurves observed from this system (Fig. 4), displayed by black dashed lines.Clearly, the light curve variability is (at least in a Newtonian system) dominated by the variability of the angle-averaged temperature profile.However, the properties of observing in the Wien tail mean that the inferred temperature profile is shifted slightly to higher temperatures, as hotter disc regions contribute exponentially more to the light curves.
It is important to note that while we have generated different luminosity distributions (in Fig. 4) by varying the observing energy E l , there is of course an exact symmetry here with varying the mean disc temperature µT (i.e., by changing the mass accretion rate, or black hole mass) at fixed observing energy.What is important (as stressed in Mummery & Balbus 2022) is the ratio between the average peak disc temperature and the observing energy kµT /E l .A ∼ 1% temperature variability can, if observed deep enough into the Wien-tail, result in order of magnitude luminosity variability.

The effects of changing the disc turbulence physics
As the analytical model of Mummery & Balbus (2022) reproduces well the full numerical stochastic-viscosity hydrodynamics models of Turner & Reynolds (2023), we now use this theory to examine the effects different assumptions about the properties of the disc turbulence have on the observed lightcurves of accretion flows at high photon energies.
There are three principal components of the viscous parameterisation which can be modified: the driving timescale 1/ω; the amplitude of variability ⟨β 2 ⟩; and the coherence length of the variability lc.Two of these, the timescale and amplitude of variability, result in relatively trivial (though important) modifications to the results.Namely, modifying the amplitude of variability ⟨β 2 ⟩ simply scales the overall temperature (and therefore luminosity) variance.This is particularly important in the Wien-tail, owing to the exponential amplification of temperature fluctuations.Similarly, modifying the driving timescale 1/ω simply scales the timescale of the peaks and troughs of the lightcurve flaring (e.g., Fig. 3), without substantially modifying the amplitude of variability.As the canonical values chosen for these two parameters are selected with reference to MHD simulations of discs (Hogg & Reynolds 2016, see further discusion in Turner & Reynolds 2023), we do not modify them further in this analysis.
The effects of modifying the coherence length of the stochastic viscosity parameterisation lc is less obvious, and physically interesting to examine.In real discs, the spatial coherence length scale lc is likely set (to first order) by the height H of the disc.The idea being that turbulent eddies will have a size given approximately by the disc thickness, and certainly will not substantially exceed this value.Turner & Reynolds (2023) presented simulations which varied H (but kept lc = H) and also ones which kept H = 0.1 but varied lc independently.For the full details on these simulations and how the coherence length is modified see Turner & Reynolds (2023) (particularly Appendix A)3 .They found that increasing the coherence length (whether or not this was accompanied by a change in H) increased the overall variability of their disc simulations.Therefore thicker discs (which are those which we would expect to be more coherent) will likely be more variable.
In order to see how this result feeds through into the high-photonenergy variability, we now consider three simulations from Turner & Reynolds (2023).Along with the fiducial simulation we have used up until now (with lc/r = H = 0.1) we consider two more which keep H = 0.1 but have lc/r = 0.15 and 0.2 respectively.While we are motivated here by considering the effect of disc thickness on variability, we have here chosen simulations with the same H.In the vertically integrated simulations of Turner & Reynolds (2023), changing H changes the inflow velocity for a given accretion rate (since a larger value of H increases the viscosity and thus lowers the surface density and increases the inflow velocity; see eq. 12).The choice of keeping H constant is therefore made to allow for easier interpretation of our results, as it ensures that the only difference between the simulations is the value of lc.Note however that, as mentioned earlier, we would expect very similar results if we considered the simulations where H was changed along with lc/r, though this is typically harder to simulate.
In Fig. 7, we show temperature snapshots for the three simulations, ordered top to bottom by increasing coherence length.Each disc temperature snapshot was made with the typical TDE parameters of M = 2 × 10 6 M⊙, Ṁ0 = 6 × 10 20 kg/s.It is relatively clear to see the effects of increasing the coherence length on the resulting disc temperature profile: unsurprisingly the disc temperature is coherent over much larger scales for longer lc (i.e. the increased coherence length in the viscosity is translated directly into that for the temperature), and the fluctuations resemble more macroscopic sheared-out spiral structures, rather than very short length scale fluctuations.
Quantitatively, increasing the coherence length both shifts the peak temperature distribution to higher temperatures, but also increases the variance in this distribution.This can be seen clearly in Fig. 8, where we display the fixed azimuth maximum temperatures Tmax, and the angle averaged maximum temperature ⟨Tmax(t)⟩ ϕ distributions, for both lc/r = 0.1 and lc/r = 0.2.
The physical reason for this change is spelled out in detail in Turner & Reynolds (2023); a disc with larger coherence length has fewer "independent" disc regions in the inner edge, with of order N ∼ 2πrI /lc independent disc regions near to the ISCO.As small scale disc structure is more highly variable (e.g., Figs. 1, 6), the smaller the number of independent regions averaged over, the larger the variance in the temperature distribution.
This modified variability coherence scale has a pronounced effect on the Wien-tail lightcurves observed from the system.In Fig. 9 we show the lightcurves for two discs of different variability coherence lengths lc/r = 0.1 (red) and 0.2 (blue).The discs and black holes of the two simulations were otherwise identical (with the same parameters as Fig. 4).Both the mean and variance of the disc luminosity are modified by changing the coherence length of the disc variability, with both increasing with lc.If we identify lc with the height of the disc, this result suggests that thicker discs will be both brighter and more variable.This may have important observational implications.
Increasing the coherence length results in both a higher mean temperature µT and temperature variance σT being inferred from fits of the analytical distribution (eq.24) to the lightcurves.In Fig. 10 we display analogous figures to Figs. 4 and 5, but now for the lightcurves generated from the lc/r = 0.2 simulation (blue points Fig. 9).All other parameters of the system are unchanged.We note a few results: the mean temperature has increased slightly (from kµT ≃ 65 eV to kµT ≃ 69 eV) which results in the larger average X-ray luminosity of the lightcurve (Fig. 9), while the temperature variance parameter has nearly doubled (from kσT ≃ 0.99 eV to kσT ≃ 1.75 eV).Perhaps surprisingly, the inferred disc size parameter has slightly reduced (from Rp ≃ 17.8rg to Rp ≃ 15rg).These inferred parameters are consistent with the shift in the angleaveraged temperature distribution shown in Fig. 8. T max (t, φ), l c /r = 0.1 . A between the fixed azimuth and angle-averaged peak temperature profiles of discs with different coherence lengths (displayed in legend; dashed curves have lc/r = 0.1, solid curves have lc/r = 0.2; both produced with M = 10 6 M ⊙ , Ṁ0 = 6 × 10 21 kg/s).Increasing the coherence length both shifts the peak temperature distribution to higher temperatures, but also increases the variance in this distribution.This change becomes clearly observable in the light curves (e.g., Fig. 9).
Importantly, the analytical model of Mummery & Balbus ( 2022) is again able to capture the observed variability of each of the different lightcurves with a simple 3-parameter model (we again find η = 5/2 best describes the data).
It is not only disc physics which can effect the observed properties of a black hole accretion flows lightcurve however.In the following section we examine how relativistic viewing angle effects may enhance (or suppress) lightcurve variability.

RELATIVISTIC VIEWING EFFECTS
The stochastic-viscosity hydrodynamic simulations of Turner & Reynolds (2023) are performed using a Newtonian theory of gravity.Naturally however, in real black hole systems there will be important relativistic effects which will modify the properties of observed accretion disc variability.These modifications will come in two distinct forms: the local accretion disc physics will be modified by relativistic gravity (this is particularly relevant for discs tilted out of the equatorial plane, for example White & Quataert 2022), and in addition there will be relativistic viewing angle effects which also modify the emergent spectrum.These viewing effects include the Doppler-shifting of the emitted photons (owing to the relativistic orbital motion of the fluid elements in the inner disc regions), and the gravitational lensing and red-shifting of the photons over their trajectory to the observer.
Running numerically expensive full general relativistic magnetohydrodynamic simulations of accretion flows lies well beyond the scope of this paper, and we are not able to examine the effects of strong gravity on the local (disc frame) variability of the accretion flow.Relativistic viewing effects however can be probed, provided that we assume that the simulated discs are evolving around Schwarzschild black holes.While this is of course not a self consistent model (the disc evolves in Newtonian gravity and is then observed in a full relativistic treatment), this modelling is likely to elucidate the key effects of varying the observer's inclination, which as far as the authors are aware is a calculation which has not yet  4).Both the mean and variance of the disc luminosity are strongly effected by changing the coherence length of the disc variability, with both increasing with lc.Identifying lc with the height of the disc suggests that thicker discs will be both brighter and more variable.
been performed in the literature.We will demonstrate in the following section that these effects can be pronounced, and the fractional variability of the luminosity can double by simply changing the observer's inclination.
In the following sub-section we will generalise the Newtonian analysis of the previous sections to a fully relativistic treatment, before performing a detailed analysis of the resultant spectra as a function of observer inclination.

Relativistic disc spectrum
The specific flux density FE from a source, as observed by a distant observer at rest (subscript obs), is given by defintion by Here, E obs is the observed photon energy and IE(E obs ) the specific intensity, both measured at the location of the distant observer.The integral is to be taken over the domain of the disc, denoted D.
The differential element of solid angle subtended on the observer's sky by a disc area element is dΩ obs .Since IE/E 3 is a relativistic invariant (e.g., Misner et al. 1973), we may write where the photon energy ratio factor fγ is the ratio of E obs to the emitted local rest frame photon energy Eemit: Figure 11.The "observed temperature" of the accretion disc simulations fγ T , projected into the x − y plane of the disc for different inclinations (θ obs ) of the distant observer, at a given snapshot in time.The photon ray-tracing calculations assume a Schwarzschild geometry.The distant observer is orientated at the very bottom of each figure (with θ obs = 0 • being directly above the disc plane), and the disc fluid is rotating in the counter-clockwise direction.The average mass accretion rate for each disc system is different, and was tuned so that the peak observed temperature (and X-ray luminosity) was comparable in each case.
The black hole mass is M = 2 × 10 6 M ⊙ .Fluid elements on the left hand side of the black hole are observed to be moving towards the distant observer, and get a corresponding Doppler boost to their temperature, while fluid elements on the right are Doppler dimmed.Gravitational red-shifting effects all fluid elements at the same radius equally.Not displayed here are the effects of gravitational lensing, which are taken into account in the luminosity calculations.For larger inclination the observed temperature is peaked in an increasingly small region of the disc, and local variability has a larger effect on the observed luminosity.
In this final expression p ϕ and −pt are the angular momentum and energy of the emitted photon (conveniently constants of motion in the general Kerr metric), U t and U ϕ are the time and azimuthal 4velocity components of the orbiting disc fluid.For a Schwarzschild spacetime these components are equal to where rg = GM/c 2 .The spectrum observed from an accretion flow, assuming that it emits thermally in its rest frame, is then equal to where E is the photon energy in the distant observer's rest frame.
The Newtonian expression (eq.15) then simply corresponds to taking the limit fγ → 1, dΩ obs → cos(i) r dr dϕ/D 2 .We take the standard approach to computing the relativistic spectral integral (eq.38).Starting with a finely space grid of photon impact parameters (described by cartesian coordinates bx, by) in the observer's camera plane, which is inclined at an angle θ obs from the z-axis of the disc, photon trajectories are traced back until they cross the disc plane θ = π/2, at an emission location (re, ϕe).At this location the photon energy shift factor fγ can be computed (eq.36), and a value for the local temperature T (re, ϕe, t) can be read in from the simulation domain.The differential element of solid angle subtended by each photon is dΩ obs = dbxdby/D 2 .The spectral integral (eq.38) can then be computed across the observer's camera plane.For more details on the algorithm used see Appendix A of Mummery & Balbus (2020).

4.2
The "observed temperature" fγT Equation ( 38) highlights an interesting property of (thermally emitting) relativistic systems.Namely that the combined effects of gravitational and Doppler energy shifting of photons along their trajectory to the observer result in a simple shifting of the temperature of the emitting material, T → fγT (this is precisely analogous to the cosmological red-shifting of the observed temperature of the cosmic microwave background, with fγ = 1/(1 + z), for example).
Therefore, it is the "observed temperature", the product of fγ and T , which is the relevant quantity for understanding the observed properties of disc lightcurves.Unlike in the cosmological analogy however, for black hole accretion problems this "observed temperature" shows a strong observer-inclination dependence, owing to the fact that the inner disc material is orbiting at a substantial fraction of the speed of light, and pronounced Doppler energy-shifts can occur if this motion is directed into the line of sight of the observer.
By computing energy shift maps fγ(r, ϕ) of photons emitted from a disc location (re, ϕe) and ultimately observed by a distant observer inclined at an angle θ obs from the z-axis of the disc plane, it is possible to examine the two dimensional temperature profiles of the Turner & Reynolds (2023) disc models, as observed by a distant observer.These two dimensional maps are shown in Fig. 11, we remind the reader that these maps were produced assuming a Schwarzschild spacetime.
In Fig. 11 we plot two dimensional profiles of the "observed temperature" of the accretion disc simulation, defined as the product fγT , projected into the x − y plane of the disc for different inclinations (θ obs ) of the distant observer (displayed in the upper right corner of each panel).Each snapshot is taken at the same instant in time, only the viewing angle is different between panels.The Figure 12.The effect of varying the distant observer inclination on the disc spectra and 0.3−10 keV lightcurves.In blue we show the lightcurves seen by a near edge-on observer θ obs = 85 • , while in red we display the lightcurves seen from a near face-on orientation θ obs = 5 • .The mass accretion rate of the two systems are tuned so that they have roughly the same average X-ray luminosity.Owing to Doppler boosting, and the correspondingly smaller disc region contributing to the luminosity of inclined disc systems, the variability of the X-ray lightcurves is substantially increased for near edge-on orientations.Note, as can be seen in the lower zoomed-in panel, that the lightcurves observed at different orientations in a relativistic system need not be precisely correlated with one another, although they often are.distant observer is located at the very bottom of each figure (with θ obs = 0 • being the orientation directly above the disc plane), and the disc fluid is rotating in the disc plane in the counter-clockwise direction.The average mass accretion rate for each disc system is different, and was tuned so that the peak observed temperature (and X-ray luminosity) was comparable in each case.The black hole mass is M = 2×10 6 M⊙, and for the remainder of this section we use the simulation with lc/r = 0.2 (i.e. as seen in the bottom panel of Fig. 7, which shows the temperature snapshot in a Newtonian theory), to maximise the observed variability.The effects of relativity introduce a clear asymmetry into the observed temperature profiles for inclined systems θ obs ≳ 45 • .The interpretation and physical cause of this result is simple, fluid elements on the left hand side of the black hole (in this Figure ) are moving with a significant component of their velocity pointing towards the distant observer, and get a corresponding Doppler boost to their temperature fγ > 1, while fluid elements on the right are moving away from the observer and are therefore Doppler dimmed fγ < 1. Gravitational red-shifting effects all fluid elements at the same radius equally, and has the net effect of reducing the observed temperature at all radii, but does not introduce an asymmetry to the profile.
Not displayed in this Figure are the effects of gravitational lensing (i.e., the bending of light rays which pass close to the central black hole).These effects are however taken into account in the full luminosity calculations.For larger inclination the disc temperature is peaked in an increasingly small region of the disc, with an observer in the Wien-tail only "seeing" a small wedge of the disc for θ obs = 85 • , while a face-on observer would "see" contributions from the full 2π azimuthal range of the disc.

Lightcurve variability as a function of inclination
In a relativistic system the observer's inclination angle has a large effect on the observed variability of a turbulent disc system.This is highlighted in Fig. 12.In Fig. 12 we display the full lightcurves generated from the 1600 snapshots of the disc temperature recorded in the simulations of Turner & Reynolds (2023), as observed at an  inclination of θ obs = 5 • (red curve) and θ obs = 85 • (blue curve); these inclinations correspond to the upper left and lower right disc temperature snapshots of Fig. 11.The black hole mass of both systems was M = 2 × 10 6 M⊙, while the average accretion rate Ṁ0 of the two systems were tweaked so that they had roughly the same average X-ray (0.3 − 10 keV) luminosity.The accretion rate had to be reduced by a factor of ∼ 3 over the range θ obs = 5 • − 85 • , with more edge-on discs brighter in the Wien-tail for fixed accretion rate (owing to Doppler boosting, which more than counteracts the usual reduction in the luminosity due to the cos i term in eq.15; see the upper panel for more details).
For large (near edge-on) inclination angles the disc regions with the largest observed temperatures shrinks to a small region with fluid elements moving towards the observer, and only this small region of the disc contributes to the observed Wien-tail luminosity.This means that observations of the Wien-tail are sensitive to an ever smaller region of the disc, and are therefore increasingly sensitive to the local physics of the disc variability.This results in a substantially larger fractional disc variability in the case of large inclinations (Fig. 12).Physically, variations on the order of a few coherence lengths lc show substantial variability amplitudes; when the flux is made up of an average over the entire azimuthal disc structure these local fluctuations are somewhat averaged out.When only a small wedge of the disc contributes to the flux, these variations show an enhanced effect.
An interesting result to note is that, as can be seen in the lower zoomed-in panel, the lightcurves observed at different orientations in a relativistic system need not be precisely correlated with one another, although they often are somewhat correlated.This highlights the point that different disc regions (more precisely different fractions of the disc angular extent) are being "seen" by different observers.
The model of Mummery & Balbus (2022) (eq. 24) still describes well the observed distributions of the disc luminosity, even with relativistic viewing angle dependence included.This is a result of the "observed temperature" property of thermal relativistic systems, meaning that the gross functional form of the Wien-tail luminosity is unchanged by including relativistic effects.In Fig. 13 we display the numerical distributions of the observed 0.3 − 10 keV X-ray luminosities for each of the viewing angles displayed in Fig. 11.The average mass accretion rate was tweaked in each case to keep the average luminosities roughly the same.It is clear to see that the fractional luminosity variability σL/µL increases with increasing inclination (upper left panel).In the upper right, lower left and lower right panels we show fits of the analytical probability density function (eq.24) to each of the numerical distributions, for a range of observer inclinations.The fit in each case is good, and the different distributions result simply from differing inferred disc temperature parameters, chiefly the differing temperature variance σT as a function of inclination.
In Fig. 14 we display both the fractional luminosity (upper) and disc temperature (lower) variability, as a function of observer inclination.Error bars on the lower plot denote the standard deviation of the fit posteriors.The fractional temperature variability σT /µT roughly doubles as the inclination is varied form 5 • to 85 • , as does the fractional luminosity variance which peaks at σL/µL ∼ 0.45, which represents extreme variability.

DISCUSSION -TIDAL DISRUPTION EVENTS AS A PROBE OF ISCO-SCALE DISC TURBULENCE
In this paper we have shown, provided that the evolution of real accretion flows is well approximated by a stochastic-viscosity hydrodynamic theory (as in the standard theory of propagating fluctuations), that tidal disruption events are ideal targets to study stochastic variability in accretion flows, as they may be observed to show extreme variability in the standard 0.3 − 10 keV X-ray bands.The reason that tidal disruption events represent particularly interesting sources from the perspective of variability studies is multi-faceted, and we list the principle reasons below.Firstly, and importantly, tidal disruption event discs typically have peak disc temperatures at the scale of kTp ≲ 100 eV (e.g., Bright et al. 2018;Wevers et al. 2019;Saxton et al. 2021;Mummery et al. 2023), a factor at least 3 lower than the typical lower bandpass energy of X-ray telescopes E l = 300 eV.This means that tidal disruption events are observed robustly in the Wien-tail of the disc spectrum, where any variability is exponentially amplified (Mummery & Balbus 2022).As we have demonstrated in this paper, relatively minor variability in the angle-averaged disc temperature profile (Fig. 6) can translate into factor unity variability in Wien-tail light curves (Fig. 4), up to scales even as high as an order of magnitude for inclined systems (Fig. 12).
Secondly, tidal disruption event discs are expected to be somewhat more compact and relatively thicker (H/R ≳ 0.1) than the analogous accretion states of X-ray binary and AGN discs, owing to the relatively high accretion rates of post-disruption debris fallback Rees (1988).As the coherence length of typical turbulent eddies is expected to be of order the scale height of the disc lc ∼ H, this means that tidal disruption event disc turbulence may well be coherent over larger scales.Equivalently, there are fewer independent disc regions near to the ISCO (N ∼ 2πrI /lc), and the intrinsic variability of these systems is larger (Figs. 8,9,10).
In addition, the large black hole masses of typical tidal disruption event hosts M• ∼ 10 6 − 10 8 M⊙ makes the orbital timescale of the innermost disc regions t orb = r 3 /GM• ∼ O(1−10) hours.This is ideal for probing the short timescale variability with dedicated follow up X-ray observations on day-week timescales.
Finally, tidal disruption event X-ray spectra typically do not show pronounced power-law components which would otherwise mask the Wien-tail.This is in sharp contrast to both X-ray binary and AGN spectra.As it is unclear how variability of a Comptonized component tracks the intrinsic variability of the accretion flow itself, observations of the soft disc-dominated component appears the best route probing intrinsic disc turbulence.
For these reasons, we recommend dedicated X-ray observational campaigns of future tidal disruption events to aid in our understanding of disc turbulence.Conversely, this can be phrased as a prediction: we expect that future short-timescale resolved analysis of Xray bright TDEs will detect ∼ order of magnitude X-ray flares on ∼ hour-to-day timescales.
One potential complication in probing the physics of disc turbulence with tidal disruption events is the possibility that these sources may be accreting at super-Eddington rates at early times, where it is less clear to what degree the predictions of the theory of propagating fluctuations robustly hold.A simple "fall-back rate" calculation implies super-Eddington accretion in a typical tidal disruption event for the first ∼ year of its evolution (Rees 1988), and it has been suggested that at these high accretion rates viewing angle effects (chiefly whether or not the observer's line of sight passes through a thick super-Eddington torus) may dominate the observed properties of the source (Dai et al. 2018).
That being said, while disc accretion at super-Eddington rates must to some degree be more complex than the thinner discs considered in this work, the main results presented here are based on simple physical principles which should be universal, and apply more broadly than to just thin discs.It seems extremely likely that all disc accretion (super-Eddington or otherwise) will remain a turbulent process, and therefore that all disc emission profiles will be stochastic to some degree.The exponential amplification of these stochastic features by the Wien tail is dependent only on the thermal nature of this emission and the observer's bandpass energy, and not on the broader physical nature of the emitter itself.It seems likely that even at early times therefore large amplitude fluctuations should be observed in tidal disruption event X-ray light curves, while at later times these systems will approach the precise theory developed here (tidal disruption events are routinely X-ray bright for ∼ 1000 days, e.g., Bright et al. 2018).

Real tidal disruption events
Given this prediction of the theory of propagating fluctuations, it is interesting to look at the X-ray lightcurve properties of known tidal disruption events to see if any hints of the behaviour described here have been previously observed. .This source shows clear large-amplitude stochastic variability on timescales as short as the observing cadence (see inset for a zoom-in of one month of data around 230 days after first detection).The X-ray spectrum of this source is dominated by a cool kT ∼ 50 eV thermal component, as expected from the models put forward in this paper.The error bars on the data points may be smaller than the marker size.
Understanding the properties of TDE X-ray light curve fluctuations will only grow in importance as the sample size of X-ray bright TDEs increases in the coming years.Indeed, the current sample of X-ray bright TDEs already appear to show significant X-ray variability, but due to sampling issues in follow-up observations (for example, non uniform observing cadences and big gaps between observations), such behaviour is currently poorly characterised.Despite this, examples of different types of behaviour already observed include rapid (∼weeks) flaring/dipping episodes (Saxton et al. 2012;van Velzen et al. 2020), sudden (Kajava et al. 2020) or smooth (Gezari et al. 2017) X-ray brightening with significant delays with respect to the UV/optical peak.While a full analysis of the lightcurve variability of tidal disruption events is beyond the scope of this work, and is postponed for future study, we do highlight one particularly interesting source: the tidal disruption event AT2019azh.
In Fig. 15 we display the X-ray lightcurve of AT2019azh (Hinkle et al. 2021;Goodwin et al. 2022), which is an X-ray bright tidal disruption event with a cool (kT ∼ 50 eV) but thermal X-ray spectrum (Hinkle et al. 2021;Mummery et al. 2023).In the latter (t > 200 days since first detection) stages of it's X-ray lightcurve evolution AT2019azh shows extreme X-ray variability, on timescales at least as short as the observational timescale.During this phase of evolution the spectrum was very soft, and dominated entirely by a thermal component (Hinkle et al. 2021;Mummery et al. 2023).There is also a robust late time detection of an optical plateau in the AT2019azh lightcurve, indicating the presence of an accretion disc (Mummery et al. 2024).
It may well be the case that AT2019azh represents an observational "detection" of the theory of propagating fluctuations, in its simplest form.There is clearly no periodic structure in the X-ray variability of AT2019azh (ruling out, for example, a repeating tidal disruption event), which appears to be stochastic.The amplitude of variability shown by AT2019azh is extremely large, as is expected by the exponential enhancement of turbulent temperature variability.
In a follow up work we will aim to quantitatively analyse AT2019azh (and other sources like it) with the numerical and analytical models presented in this paper.

CONCLUSIONS
In this paper we have analysed the properties of accretion disc variability (driven by fluctuations in the angular momentum transport) as observed at high photon energies (i.e., photon energies larger than the peak disc temperature).We have shown that small amplitude temperature fluctuations (of order ∼ percent level) can result in up to ∼ order of magnitude luminosity variability, owing to the exponential enhancement of the Wien-tail.This variability is well described by the analytical model of Mummery & Balbus (2022).This excellent fit between analytical and numerical theory means that the (significantly cheaper to compute) analytic theory can be used to probe variability in astrophysical disc systems, and potentially elucidate the nature of disc turbulence.
By utilising the numerical stochastic-viscosity hydrodynamic models of Turner & Reynolds (2023), we have shown that both the physics of the disc turbulence (in particular the turbulent coherence length), and relativistic viewing angle effects, can have pronounced effects on the observed variability.Increasing the turbulent coherence length (likely set by the disc thickness) results in more pronounced large scale structures to develop in the flow (Fig. 7), which results in an increase in luminosity variability.Observing a relativistic disc system at a larger inclination results, at fixed X-ray luminosity, in an increase in the observed variability of the disc.This is a result of the "zooming in" onto smaller scale structure in the small disc region which is Doppler boosted towards the observer, which dominates the total flux.As far as the authors are aware this is the first such analysis of relativistic viewing effects.At fixed inclinations we have demonstrated that fainter disc systems have a larger amplitude of variability (Fig. 4; see also Mummery & Balbus 2022).
Finally, we argue that tidal disruption events may represent the ideal astrophysical systems with which to observationally probe the small length scale turbulent physics of black hole accretion.Tidal disruption event discs have the fortunate properties of (i) being typically much cooler than observing bandpasses, putting their emission in the Wien tail; (ii) likely having thicker discs (and larger turbulence coherence lengths) than soft-state X-ray binaries; (iii) lacking complicating power-law components in their X-ray spectra, allowing a cleaner look into the inner disc, and (iv) naturally evolving on ∼ hour-to-day timescales, owing to their large black hole masses.Indeed, the tidal disruption event AT2019azh (Hinkle et al. 2021;Goodwin et al. 2022) may well already have been observed to show properties predicted by the theory of propagating fluctuations.Computing Service (www.csd3.cam.ac.uk), provided by Dell EMC and Intel using Tier-2 funding from the Engineering and Physical Sciences Research Council (capital grant EP/T022159/1), and DiRAC funding from the Science and Technology Facilities Council (www.dirac.ac.uk).

Figure 2 .
Figure2.A set of 1600 observed disc spectra, generated from each of the disc temperature snapshots (see Fig.1).The black hole mass was set to M = 2 × 10 6 M ⊙ , and the average mass accretion rate Ṁ0 = 6 × 10 20 kg/s, so as to produce typical 0.3 − 10 keV X-ray luminosities typical of a tidal disruption event (a particularly relevant astrophysical system).It is clear that the fractional variability in the disc luminosity increases as a function of observing energy, and is particularly pronounced in the Wien tail.The conventional 1D disc theory spectrum for these parameters is shown by a black dashed curve.Note that for the purposes of computing light curves the disc spectrum is evaluated above 1 keV, but this is not shown here.

Figure 3 .
Figure 3. Example lightcurves observed across 0.3 − 10 keV, for the spectra depicted in Fig. 2. Upper: the full lightcurve, while a temporal zoom-in is shown in the lower panel.For a black hole mass of M = 2 × 10 6 M ⊙ , the gravitational timescale tg = GM/c 3 ≃ 10 seconds, and so each point here is separated by ∆t = 10 3 tg ∼ 2 hours 45 minutes.Clearly, when observed in the Wien tail, order of magnitude X-ray luminosity fluctuations are possible on ∼ hour timescales.
Figure 4.The numerical (coloured histograms) and theoretical (purple curves) distributions of the X-ray lightcurves.Each lightcurve was generated from the same disc (with M = 10 6 M ⊙ , Ṁ0 = 6 × 10 21 kg/s), but was observed with a different lower bandpass energy E l .The blue (brightest) histogram represents the lowest observing energy, while the red (dimmest)shows the highest.These different bandpass energies were chosen so that the average luminosity across the bands varied by a factor ∼ 100.The purple curves are simultaneous fits of the analytical distribution (eq.24) to the numerical data; the fit is excellent.

Figure 7 .
Figure 7.The effects of increasing the coherence length of the stochastic viscosity prescription on the disc temperature profile: upper panel lc/r = 0.1, middle panel lc/r = 0.15, lower panel lc/r = 0.2.As the coherence length of the viscosity is increased (from upper to lower panels) the temperature variability becomes increasingly structured on large scales, and clear spiral structures are visible.This results in larger amplitude luminosity variability, particularly when observed in the Wien-tail.

Figure 9 .
Figure9.Example lightcurves for discs of different variability coherence lengths lc/r (displayed on plot).Upper: the full lightcurve, while a temporal zoom-in is shown in the lower panel.The discs and black holes of the two simulations were otherwise identical (with the same parameters as Fig.4).Both the mean and variance of the disc luminosity are strongly effected by changing the coherence length of the disc variability, with both increasing with lc.Identifying lc with the height of the disc suggests that thicker discs will be both brighter and more variable.

Figure 10 .
Figure10.Upper: Numerical (coloured histograms) and theoretical (purple curves) distributions of the X-ray lightcurves.Each lightcurve was generated with the same disc parameters as in Fig.4(i.e., M = 10 6 M ⊙ , Ṁ0 = 6 × 10 21 kg/s), except for the coherence length of the stochastic viscosity, which was doubled.The coherence length is anticipated to be of order the scale height of the disc, and so this represents a thicker disc model.The observing bands were kept identical to Fig.4.The blue (brightest) histogram represents the lowest observing energy, while the red (dimmest) shows the highest.The purple curves are fits of the analytical distribution (eq.24) to the data; the fit is excellent.Lower: correlations between the parameters of the Xray luminosity probability density function, and their posterior distributions.Doubling the coherence length of the stochastic viscosity has approximately doubled the fractional temperature variability σ T /µ T ∼ 2.5%, and resulted in significantly enhanced luminosity variability.

Figure 13 .
Figure13.The X-ray luminosity distributions of relativistic accretion flows observed at increasing inclinations.Upper left: the numerical distributions of the X-ray luminosity for inclinations in the range θ obs = 5 • − 85 • .The mass accretion rate of the different inclinations were tuned to give roughly the same average X-ray luminosity.The increasing luminosity variance with inclination is clear to see.Upper right, lower left, and lower right: individual distributions for three different inclinations, and fits to the observations of the analytical distribution (eq.24; purple curves).The analytical results well describe the distributions, even with relativistic viewing effects included.The observed temperature variance is a function of inclination, as shown in Fig.14.

Figure 14 .
Figure14.Upper: the fractional X-ray luminosity variability as a function of observer inclination.Lower: the inferred fractional temperature variability, from fits of the probability distribution (eq.24) to the lightcurves, also as a function of observer inclination.Error bars denote the standard deviation of the fit posteriors.The fractional temperature variability roughly doubles with inclination, as does the fractional luminosity variance.

Figure 15 .
Figure15.The X-ray (0.3 − keV) lightcurve of the tidal disruption event AT2019azh.This source shows clear large-amplitude stochastic variability on timescales as short as the observing cadence (see inset for a zoom-in of one month of data around 230 days after first detection).The X-ray spectrum of this source is dominated by a cool kT ∼ 50 eV thermal component, as expected from the models put forward in this paper.The error bars on the data points may be smaller than the marker size.

ACKNOWLEDGMENTS
AM was supported by a Leverhulme Trust International.SGDT acknowledges support under STFC Grant ST/X001113/1, and previously under an STFC PhD Studentship.For the purpose of Open Access, the authors have applied a CC BY public copyright licence to any Author Accepted Manuscript version arising from this submission.The hydrodynamical simulations used in this work (which were originally presented in Turner & Reynolds 2023), were performed using resources provided by the Cambridge Service for Data Driven Discovery (CSD3) operated by the University of Cambridge Research