Truncated, Tilted Discs as a Possible Source of Quasi-Periodic Oscillations

Many accreting black holes and neutron stars exhibit rapid variability in their X-ray light curves, termed quasi-periodic oscillations (QPOs). The most commonly observed type is the low-frequency ( ≲ 10 Hz), type-C QPO, while only a handful of sources exhibit high frequency QPOs ( ≳ 60 Hz). The leading model for the type-C QPO is Lense-Thirring precession of a hot, geometrically thick accretion flow that is misaligned with the black hole’s spin axis. However, existing versions of this model have not taken into account the effects of a surrounding, geometrically thin disc on the precessing, inner, geometrically thick flow. In Bollimpalli et. al 2023, using a set of GRMHD simulations of tilted, truncated accretion discs, we confirmed that the outer thin disc slows down the precession rate of the precessing torus, which has direct observational implications for type-C QPOs. In this paper, we provide a detailed analysis of those simulations and compare them with an aligned truncated disc simulation. We find that the misalignment of the disc excites additional variability in the inner hot flow, which is absent in the comparable aligned-disc simulations. This suggests that the misalignment may be a crucial requirement for producing QPOs. We attribute this variability to global vertical oscillations of the inner torus at epicyclic frequencies corresponding to the transition radius. This explanation is consistent with current observations of higher frequency QPOs in black hole X-ray binary systems.


INTRODUCTION
Quasi-periodic oscillations (QPOs) are among the most intriguing phenomena observed in the X-ray light curves of accreting black hole and neutron star X-ray binary systems.In black holes these features are identified by relatively broad peaks in the power spectrum, where the broadening may result from the modulation of the oscillation frequency itself, or the finite lifetime of the oscillation.The frequency of these QPOs is believed to be closely linked to the strong gravity of the compact object (van der Klis 2000; Kluzniak & Abramowicz 2001b;van der Klis 2004;Ingram & Motta 2019).Therefore, understanding the physical nature of QPOs can be a step towards a better understanding of compact objects.Furthermore, timing studies probe smaller physical scales (of order the size of the compact object) than is possible with direct imaging of black hole binaries-only in the case of nearby supermassive black holes can such scales be probed with the current generation of the Event Horizon Telescope (EHT).
In addition to the rapid variability, a fair population of neutron star and black hole X-ray binaries also exhibit "state transitions", a phenomenon in which the spectral state of the source evolves concurrently with its luminosity, often tracing a "q"-shape in a hardnessintensity diagram (Narayan et al. 1996;Esin et al. 1997;Belloni 2005;McClintock & Remillard 2006;Begelman & Armitage 2014).The two major spectral states observed in these systems, along with ★ Contact e-mail: deepika@mpa-garching.mpg.de the various intermediate states, are the "soft" and "hard" states.In the soft state, the X-ray spectra are predominantly composed of a soft, thermal component (below approximately 3 keV).This soft emission likely arises from a relatively cold, optically thick and geometrically thin disc, exhibiting thermal blackbody-like radiation.On the other hand, in the hard state, the spectra are dominated by a power-law component (extending roughly between 10 -100 keV).This powerlaw emission results from the Compton upscattering of soft seed photons emitted from the disc by a radiatively inefficient cloud of hot electrons, often referred to as the "corona."The location and geometry of the corona are still the subject of considerable debate; some models suggest that it resides at the base of the jet (Martocchia & Matt 1996;Fabian et al. 2009;Parker et al. 2015), while others consider it to sandwich the disc (Galeev et al. 1979;Haardt & Maraschi 1993;Svensson 1996;Beloborodov 1999).Alternatively, the truncated-disc model envisions a third possible geometry with a standard-thin disc truncated well outside the last stable orbit, and the corona as a geometrically thick, radiatively inefficient accretion flow filling the inner gap (Eardley et al. 1975;Esin et al. 1997;Liu et al. 2007).Despite the uncertainty regarding the location of the corona, all of these models can reproduce the observed X-ray spectra during the "state transitions" by associating them with appropriate changes in the geometry and strength of the accretion disc and corona, as well as the presence or absence of other physical components such as a jet.
QPOs of various flavours are observed during the different spec-tral states of a given source.Based on the frequencies, they are broadly classified into high-frequency (≳ 60 Hz) and low-frequency (≲ 30 Hz) QPOs (Remillard & McClintock 2006;Belloni 2010).High-frequency QPOs (HFQPOs) are typically observed in the highflux/luminous accretion states (Belloni et al. 2012;Motta 2016).If occurring in the thin disc, they may be associated with trapped gravity (-) and corrugation (-) modes within the inner regions of the disc (Kato & Fukue 1980;Wagoner 1999;Silbergleit et al. 2001), or the accretion disc responding quasi-periodically to a wobbling jet through the propagation of fast-magnetosonic waves (Ferreira et al. 2022), or the oscillations of spiral density waves moving between the disc's inner radius and the Lindblad radius at some frequency (Tagger & Pellat 1999;Rodriguez et al. 2002).If they occur in the geometrically thick (toroidal) flow, the HFQPOs may correspond to modes occurring at combinations of the orbital and epicyclic frequencies at characteristic radii (Rezzolla et al. 2003;Bursa et al. 2004;Blaes et al. 2006;Fragile et al. 2016) or other gravity/inertial pressure modes within the accretion disc (Abramowicz et al. 2006).
Low-frequency QPOs (LFQPOs) are further classified into type-A, B & C, of which the first two are observed in the intermediate states during the transition to the soft spectral state (Motta 2016).Type-C QPOs are commonly observed in most accretion states (Motta et al. 2012;Motta 2016), although they are particularly prominent in the hard spectral state.For these QPOs, both the fractional QPO amplitude and the phase-lag between the hard and soft photons, measured at the QPO frequency, exhibit strong dependence on the inclination angle (Motta et al. 2015;van den Eĳnden et al. 2017).This inclination dependence provides compelling evidence for a geometric origin for this QPO.For a more comprehensive review on type-C QPOs, we recommend Ingram & Motta (2019).Most of the models invoked to explain LFQPOs are based either on the geometry or the intrinsic properties of the accretion disc.Several potential mechanisms have been suggested, including oscillations of a standing shock in the accretion flow have been proposed as a source of LFQPOs.However, for such a shock to generate the required frequencies, it would need to be located at a considerable distance from the black hole (Chakrabarti & Titarchuk 1995).
Another popular model is the relativistic precession model, which associates observed QPOs with different frequencies derived from geodesic particle motion in general relativity, notably, nodal and periastron precession, as well as the particle orbit itself (Stella & Vietri 1998;Stella et al. 1999).A variant of this model, which appears to align well with the observations of the type-C QPO, entails the precession of a finite inner region of the accretion flow instead of an infinitesimal particle or ring (Ingram et al. 2009).This model corresponds with the truncated disc framework of the hard accretion state, as the precession is attributed to the hot, thick flow within the truncation radius (i.e., the corona) (Ingram & Done 2011).This model offers several advantages, including lowering the precession frequency into the appropriate range, naturally explaining the rise in QPO frequency during the early phases of an outburst cycle (when the truncation radius is moving inwards), and readily reproducing the inclination dependence observed in QPO studies (since precession is more apparent the closer to the orbital plane one observes).Furthermore, this model inherently accounts for the observed correlation between the QPO phase and the centroid of the iron emission line, when studied with phase-resolved spectroscopy (Ingram et al. 2016).
Recently, the rigid-body precession model of the torus has been extended to also explain the HFQPOs with a 3:2 ratio that are occasionally observed in the same sources (Fragile et al. 2016).In this extended framework, the HFQPOs are manifestations of natural, global oscillation modes of the same hot, thick flow that is under-going Lense-Thirring precession.For both the HFQPOs and type-C QPO, the location of the truncation radius plays a critical role in setting their frequencies.During an X-ray outburst, as the source luminosity increases, the truncation radius moves inward, causing all QPO frequencies to rise.Eventually, the QPOs vanish when the truncation radius reaches the last stable orbit in the soft state.However, it is essential to note that the key assumption of this QPO model, particularly in relation to the type-C QPO, is that the inner, hot, and geometrically thick accretion flow can precess as a rigid body independently of the surrounding thin disc.
Previous numerical work employing magnetohydrodynamic simulations of isolated, tilted thick discs have demonstrated global Lense-Thirring precession at frequencies consistent with observations (Fragile et al. 2007;Liska et al. 2018;White et al. 2019).Additionally, independent analytical and numerical studies have revealed that isolated gas tori exhibit global oscillation modes that match the general properties of observed HFQPOs (Blaes et al. 2006;Mishra et al. 2017).However, none of these studies have taken into account the full truncated-disc geometry, especially the impact of the surrounding geometrically thin disc on the QPO frequencies.A particular concern is understanding how the outer, thin disc affects the precession of the inner torus.To address this, we conducted a comprehensive set of 3D general relativistic magnetohydrodynamics (GRMHD) simulations of tilted and truncated discs around black holes using Cosmos++ (Anninos et al. 2005).Through these simulations, we demonstrated a significant decrease in the precession rate of the inner torus due to the exchange of angular momentum between the inner torus and the outer thin disc (Bollimpalli et al. 2023).This finding bears substantial implications for studies utilizing QPO frequencies to infer the truncation radius, as well as the mass and spin of the black hole.The purpose of this paper is to delve into further details regarding these simulations, particularly focusing on the properties of the accretion flow and the presence of additional QPOlike features.In addition to the simulations of tilted, truncated flows from Bollimpalli et al. (2023), we include a simulation of an aligned, truncated flow around a black hole to facilitate direct comparison.
The remainder of this paper is organized as follows: In Section 2, we describe the initial setup of our simulations, including both aligned and misaligned configurations with respect to the black hole's spin axis.We then delve into the results of our simulations, particularly focusing on the properties and dynamics of the accretion disc (Section 3).In Section 4 we present the temporal variability observed in our simulations.Subsequently, in Section 5, we provide a contextual framework for our findings by comparing them to previous works and establishing their relevance to observations.Finally, we provide our concluding remarks in Section 6.

Physical setup
We conducted four, fully 3D GRMHD simulations of truncated accretion discs, as listed in Table 1: Two high-resolution, 4-level simulations, aimed at studying tilted and untilted scenarios, and two 3-level, low-resolution simulations to explore different black hole spins.The two 4-level simulations, a9b0L4 and a9b15L4, differed in terms of the initial alignment of angular momentum axes between the disc and the black hole (with dimensionless spin parameter  * = / = 0.9 in both cases).In the former, the axes were aligned, while in the latter, they were misaligned by an angle  0 = 15 • .In the tilted disc simulations, the initial disc angular momentum vector was oriented Figure 1.Initial setup.Top: Pseudocolour plot of gas density at  = 0, overlaid with our 4-level mesh.Bottom: Pseudocolour plot of the ratio of magnetic pressure to gas pressure overlaid with magnetic field contours.The units of the axes are in  / 2 , while the density is plotted in code units.
along the -axis, while the black hole spin axis was tilted in the −direction.The 3-level simulations, a5b15L3 and a9b15L3, were both tilted 15 • and varied only in their black hole spins, with  * = 0.5 and 0.9, respectively.
The initial setup of all simulations comprises a finite torus surrounded by a thin disc, as depicted in Figure 1.The torus is initialized using the method described by Chakrabarti (1985), assuming a constant specific angular momentum.An adiabatic equation of state is assumed with a polytropic index of  = 5/3.We chose the inner edge of the torus at  in = 6.5  / 2 and pressure maximum at  cen = 9  / 2 .The torus is surrounded by a thin slab with a fixed height of  = 0.4  / 2 , extending from 15  / 2 to the outer boundary of the simulation domain at approximately 250  / 2 .The density distribution within the slab is determined by where  =  sin  is the cylindrical radius and  g =  / 2 .The angular momentum within the disc follows a Keplerian distribution.To preserve the desired thin structure of the outer disc, we implement an artificial cooling function (only in the thin disc region), as outlined in Fragile et al. (2012a), which is applied exclusively to radii beyond 15  / 2 .
Both the torus and the slab are threaded with numerous small poloidal loops of magnetic field with alternating polarity, starting from a magnetic vector potential of the form: where  g represents the local gas pressure and  = 10[(/) 2 +  2 /( − ms ) 2 +  2 /(40 1.5 − ) 2 −1].Here  = 0.6 and 0.4  / 2 in the torus and disc regions, respectively, while  ms signifies the location of the marginally stable orbit.The resulting magnetic field is normalized such that the initial ratio of the gas pressure to the magnetic pressure is ⩾ 10 throughout both the torus and disc, as shown in Fig. 1.This magnetic field configuration impedes the accumulation of strong net flux in the inner regions, thus ruling out the formation of a magnetically arrested disc (MAD), which is appropriate for the  (Fragile et al. 2023).The background region surrounding the disc and torus is initialized with a non-magnetic gas with low-density,  = 10 −5  max  −1.5 , and internal energy density,  = 10 −7  max  −2.5 .Here  max and  max correspond to the maximum gas density and internal energy density within the disc/torus.
The simulations are performed in modified Kerr-Schild sphericalpolar coordinates, enabling the placement of the inner domain boundary inside the black hole event horizon.The grid is logarithmically spaced in , covering the range  ∈ [1.4  g , 40 1.5  g ], and uniformly spaced in , encompassing the entire [0, 2] domain.The polar angle  is uniformly spaced in the coordinate  2 , with the grid spacing defined as where ℎ = 0.5 is used to concentrate the grid cells close to the midplane (McKinney 2006).To avoid computing metric terms at the poles, a small cone of opening angle 10 −15  is excised.
All simulations have a base resolution of 48 × 32 × 32, accompanied by two or three additional levels of static mesh refinement.Adjacent refinement levels differ in resolution by a factor of 2 in each dimension, resulting in a fiducial resolution of either 192 × 128 × 128 or 384 × 256 × 256 for most of the grid.The initial two levels of refinement are utilized to adequately resolve both the torus and the thin slab, which is common for both the 3-level and 4-level simulations.In the 4-level simulations, an additional refinement layer is implemented to enhance the resolution of the thin disc over the range  ∈ [10.5,  g , 60,  g ].In the untilted disc simulation, the region within 1.4 < / g < 21 and 0.28 <  < 0.72 undergoes two levels of refinement.However, for the tilted disc simulation, the second level of refinement for the torus region is adjusted to 1.4 < / g < 24 and 0.2 <  < 0.8.The thin disc region is refined by incorporating two levels of refinement within 12.6 < / g < 253 and 0.45 <  < 0.55, with the possibility of a third layer, as described earlier.
Outflow boundary conditions 1 are applied at both the inner and outer radial boundaries.Near the poles, pole axis boundary conditions are utilized, where information from the corresponding zone across the pole is used for calculating gradients, although fluxes, magnetic fields, and electric fields (emf's) are set to zero on the faces and edges that touch the pole.Additionally, periodic boundary conditions are implemented along the azimuthal direction.To maintain the zero-divergence of the magnetic fields, the vector potential method (Fragile et al. 2019) is employed in the untilted disc simulation, while in the tilted disc simulations, the constrained transport scheme (Fragile et al. 2012b) is utilized.Whenever density or pressure is boosted as a result of floors and ceilings, the boost is performed in the drift frame to preserve momentum along the magnetic field lines (Ressler et al. 2017).
As expected, once the simulations are initialized and evolve, the initial magnetic fields trigger the magnetorotational instability, leading to turbulence that facilitates the transport of angular momentum and enables accretion.The energy generated in this process causes the torus to inflate, resulting in the formation of a hot and thick flow in the inner regions.The cooling function employed in the thin slab region ( ⩾ 15  / 2 ) acts as a heat sink and helps to regulate the thickness of the disc at the target value of / ≈ 0.05.The simulations have been run for over 25, 000  / 3 for the 4level simulations and approximately 15, 000  / 3 for the 3-level simulations.In all cases, the simulations successfully evolve into a two-component flow, as anticipated.

RESULTS
At a qualitative level, the results of both the 4-level and 3-level simulations exhibit similarities, which is promising as it indicates that our final conclusions are independent of the resolution.A significant finding from this study is that the simulations of tilted and untilted discs demonstrate notable and intriguing differences in both their dynamics and variability characteristics.These distinctions will be thoroughly explored in the subsequent two sections.

Accretion flow properties
We initially aim to confirm that the simulations exhibit the desired two-component structure.One approach is to examine the densitysquared weighted scale height of the disc, computed as: where  denotes the fluid density,  is the metric determinant, and  mid represents the angular position of the disc's mid-plane.In the case of an untilted disc,  mid is taken to be /2.However, in the case of a precessing disc, since the disc mid-plane departs from the symmetry plane of the grid, we need an alternate procedure.We choose to transform disc mid-plane back to the - coordinate plane using the coordinate transformation described in White et al. (2019).In the remainder of this paper, whenever used, we denote the quantities in these transformed coordinates by a prime symbol ( ′ ).Fig. 2 illustrates the evolution of the disc scale height as a function of radius and time for both the untilted and tilted disc simulations.In both cases, a distinct two-component flow is evident, comprising a geometrically thick disc in the inner region, surrounded by a geometrically thinner disc.The geometrically thin disc beyond 15  / 2 maintains a scale height close to 0.05 throughout the entire simulation, confirming that our cooling function in the thin disc region (> 15  / 2 ) operates as expected.Also, the inner thick flow, which is at least three times thicker than the outer thin disc, maintains its scale height due to the advective nature of the thick accretion flows; most of the dissipated energy is advected inwards along with the matter (and eventually onto the black hole) before altering the vertical structure of the flow locally.Note that, in the tilted disc simulations, the scale height of the torus region is noticeably thicker compared to the untilted disc simulation.This increase in scale height is likely a consequence of additional dissipation associated with tilted discs (Dexter & Fragile 2013), which agrees with our observation of at least an order of magnitude higher temperatures in the tilted simulation (a9b15L4) compared to the untilted one (a9b0L4).
In the top panel of Fig. 3, we show the time evolution of the mass accretion rate, measured in arbitrary units, for the high-resolution tilted simulation, a9b15L4.We define this quantity as such that a negative sign implies an inflow of matter toward the black hole.Here   is the radial component of the fluid 4-velocity,   .In general, the accretion rate in the thick disc region is at least an order of magnitude lower than that in the outer, geometrically thin region.The mismatch in  implies that material must be accumulating near the transition region (roughly between 11 and 19  / 2 ), which is consistent with our findings of increased surface density at those radii (see for e.g., Fig B1).The mass accretion rate gradually increases in both the thick and thin-disc regions throughout the simulation.Similar results are observed in our other tilted disc simulations, as well.By closely comparing this figure to the bottom panel of Fig. 2, one can find correlation between the changes in the accretion rate and the scale height.
The bottom panel of Fig. 3 shows the mass accretion rate as a function of time at  = 5 and 25  / 2 for all four simulations.Upon comparing the accretion rate at 5  / 2 for simulations a9b15L3 and a5b15L3, we find that it is higher for the higher spin case.This could be due to the enhanced angular momentum transport caused by standing-shock features, which become more prominent for higher spins.These standing shocks are characteristic to tilted, thick accretion flows, and are formed along the line of nodes between ), scaled to the local orbital period, for all three tilted disc simulations.The Lense-Thirring (LT) precession timescales are shorter than the accretion timescales in the thick disc region (6 ≲ / g ≲ 20), while the sound-crossing timescales are shorter than the precession timescales everywhere except in the very inner regions of the thick discs, thus allowing for rigid-body precession of the thick disc.
the black hole symmetry plane and the disc midplane (Fragile & Blaes 2008;Generozov et al. 2014;White et al. 2019).More notable is that the accretion rates for the misaligned discs exhibit oscillations that stand out clearly in the  = 25  / 2 panel of Fig. 3, while such oscillations are absent for the untilted simulation.We will elaborate on this in Section 4, when we discuss time variability.
An important diagnostic for tilted accretion discs is comparing characteristic timescales.Fig. 4 illustrates the accretion,  acc = /⟨  ⟩, sound-crossing,  cs = 2/⟨  ⟩, and Lense-Thirring precession timescales,  LT = 2/Ω LT , all normalized by the dynamical timescale,  dyn = 2/⟨Ω⟩, for all three tilted simulations.Here, ⟨  ⟩, ⟨  ⟩ and ⟨Ω⟩ are the density-weighted averages of radial velocity, local sound speed and the angular velocity, respectively and the time averaging is done from  = 10000  / 3 till the end of the simulation.Also, is the Lense-Thirring precession rate, with Ω the orbital angular velocity, and Ω ⊥ the vertical epicyclic one.For precession to occur, we need the precession timescale to be smaller than the accretion one (Bardeen & Petterson 1975).Otherwise, the accreted angular momentum can impede precession and lead to disc alignment.Fig. 4 reveals that this criterion is only fulfilled in the thick disc region (6 ≲ / g ≲ 20), indicating precession should only be expected in that region.Furthermore, we find that the sound-crossing timescale is shorter than the precession timescale at nearly all radii (excluding the innermost regions), causing a strong coupling of the disc material.This leads us to our anticipation of rigid-body precession (Larwood et al. 1996).We present more details on this in the following section.

Tilt & Twist
To investigate the precession and warp in tilted disc simulations, it is useful to evaluate the Euler angles,  and , which represent the tilt and twist, respectively.Following the approach introduced in Fragile et al. ( 2007), we define the tilt angle as the measure of misalignment between the disc and the black hole spin axis as where is the angular momentum vector of the black hole, and for  = 1, 2, 3, which correspond to the Cartesian vector components of the disc's angular momentum vector measured in asymptotically flat space.Here, and Here   =  h + 2 m     +  g +  m   −     is the stressenergy tensor of the fluid in the magneto-hydrodynamic limit, where h = 1 +  +  g / is the specific enthalpy,  is the specific internal energy density,  m is the magnetic pressure, and   is a four-vector magnetic field measured by an observer in the frame co-moving with the fluid, and relates to the divergence-free spatial magnetic field vector through   = √ −(   0 −  0   ).At  = 0, the angular momentum vector of the disc is aligned along the direction of the ẑ unit vector.Fig. 5 shows the evolution of the tilt angle,  (measured in degrees), for each of the tilted simulations at three different radii: 5, 10, and 25  / 2 .In all simulations, we find that the thick-disc region closest to the black hole actually initially tilts further away from the black hole before eventually stabilizing at around 35 • (see the top panel in Fig. 5).This tendency for tilted thick discs to tilt away from the black hole symmetry plane at small radii is a common observation in numerical simulations (Fragile et al. 2007;Liska et al. 2018;White et al. 2019) and is consistent with expectations in the bending-wave limit (Ivanov & Illarionov 1997).At slightly larger radii (middle panel of Fig. 5), the evolution of the tilt angle exhibits some differences between  * = 0.5 and 0.9 simulations.While  initially increases in all simulations, it quickly reaches a stationary value for both the lowand high-resolution  * = 0.9 simulations.In contrast, the tilt in the  * = 0.5 simulation appears to continue growing over time.We  believe that this initial surge in  is rather a consequence of the torus responding to the bending waves propagating within it (Zhuravlev & Ivanov 2011;Zhuravlev et al. 2014).
The decrease in the tilt angle with increasing radius highlights the warped nature of these discs.This can be visualized in Fig. 6, where the bottom panel illustrates the time evolution of tilt as a function of radius for simulation a9b15L4.The top panel displays a radial profile of the same quantity, time-averaged over the interval 10000 − 25000  /3 .In the inner regions, within 20  /2 , the tilt angle increases as the radius decreases before subsequently decreasing very close to the black hole.This behaviour remains relatively stable over time.The radial profile of ⟨⟩ t reveals that the thin disc beyond 20  / 2 maintains its initial misalignment at  0 = 15 • , while the inner thick disc and transition region exhibit the expected bending wave pattern (Ivanov & Illarionov 1997).
We now calculate the precession angle, , which measures the amount of twisting or precession of the disc's angular momentum vector around the black hole spin axis.The precession angle is calculated as We also consider the projection of J BH × J disc () onto x to avoid the degeneracy in cosine for angles greater than 180 • .Fig. 7 illustrates the precession rate, /, for the simulation a9b15L4.As previously shown in Fig. 1 of Bollimpalli et al. (2023), a precession front originating from the inner thick-disc region propagates outward, causing the disc to warp as the front expands over the bending wave timescale 2 (depicted by the black dashed curve).The uniform colour in the precession rate within the 6 − 15  / 2 range suggests rigid-body precession 3 of the thick-disc region.It is evident from this figure that the precession rate of the inner region decreases significantly when the bending wave reaches the transition region (close to 15  / 2 ).As we previously reported in Bollimpalli et al. (2023), the angular momentum fluxes from the outer thin disc, counteract the Lense-Thirring torque that drives the precession of the inner torus, resulting in a deceleration of its precession.
As the bending wave passes through the thin disc, it also undergoes a brief period of precession.This happens for two reasons: a) the accretion timescales are longer than the precession timescales in the disc before the bending wave reaches it, and b) in all our simulations, the effective stress parameter 4 , , remains less than the disc scale height everywhere.The latter could be a consequence of underresolving the outer thin disc, although conducting higher-resolution simulations to test this hypothesis would be computationally expensive.However, once the bending wave has propagated through the thin disc, the accretion timescales become shorter than the precession timescales.Consequently, the thin disc eventually ceases to precess, and any warps are diffused by dynamical processes within the disc.Our simulation of an isolated torus with the same initial setup as the torus in simulation a9b15L4 exhibits the expected precession rate for an isolated torus (see Appendix A for further details).

TIME VARIABILITY
In this section, we perform a comprehensive timing analysis on various accretion flow quantities to investigate the presence of QPOs.To compute the power spectra, we employ the Lomb-Scargle periodogram method (Lomb 1976;Scargle 1982;VanderPlas 2018), focusing on density-weighted averages of a quantity, , calculated as follows: where (, , , ) represents the density distribution.The power spectra are calculated within a frequency range bounded by the minimum frequency set by the inverse of the time duration being analyzed and by the maximum frequency of 1/(2Δ), where Δ is the time interval between successive data dumps, which is provided in the last column of Table 1.Fig. 8 displays the radially-dependent power spectra of the densityweighted average of the poloidal component of the transport velocity 5 ,  , over a spherical shell for all simulations.In the case of the untilted disc (upper-left panel), there is clear power along the epicyclic frequency in the thin disc region; however, there is no evident global QPO-like feature, which would appear as a horizontal band of power.
In contrast to the untilted disc simulation, the power spectral densities (PSDs) in the tilted simulations exhibit clear signatures of variability at at least two discrete frequencies each, ranging between 0.001 and 0.006  3 / , with all three simulations exhibiting oscillations at a frequency close to 0.002  3 / .Given that we are looking at power in the vertical component of velocity, these features suggest that the inner disc undergoes vertical oscillations.We will provide a detailed discussion on the nature of these oscillations in the remainder of this section.
It is worth noting the enhanced variability power near the vertical epicyclic frequency curve (Ω z , dotted black) or, more accurately, between the vertical and radial epicyclic frequency (, dashed black) curves.Such variability at local epicyclic frequencies is expected in geometrically thin discs (Reynolds & Miller 2009;Mishra et al. 2020), but also in tilted thick discs, as they exhibit eccentric orbital trajectories (Ferreira & Ogilvie 2008;Henisey et al. 2009).
moving frame to the total pressure.Both the stress and the pressure terms include contributions from the magnetic field. 5The transport velocity,   , is related to the fluid-four velocity through

Trapped 𝑔-modes
Trapped -modes6 in the inner regions of accretion disc have been studied as a possible source of high-frequency QPOs for a long time (Kato & Fukue 1980;Okazaki et al. 1987;Nowak & Wagoner 1992;Ferreira & Ogilvie 2008).While these modes are observed in hydrodynamic simulations (Reynolds & Miller 2009;Dewberry et al. 2020a;Mishra et al. 2020), they are not commonly observed in magnetohydrodynamic (MHD) discs, possibly due to damping effects associated with the magnetorotational instability (MRI).However, recent studies by Dewberry et al. (2020b) have shown that sufficiently strong eccentric distortions can excite trapped inertial waves in MHD simulations.Tilted discs, which exhibit eccentric orbits driven by the tilt (Fragile & Blaes 2008), offer a favourable environment for the survival of these trapped g-modes even in the presence of MRI turbulence.
In Fig. 9, we present the dynamical power spectrum of the same quantity, ⟨  ⟩  ,  , computed at  = 2  / 2 .While this radius technically falls inside the ISCO for both the 0.5 and 0.9 spin cases, we deliberately chose it to effectively showcase all the QPO features across all three simulations simultaneously.To compute this, we apply a periodogram analysis to a sliding time window of 3000  / 3 .The long-dashed, black line in Fig. 9 represents the maximum of the radial epicyclic frequency ( max ), which is determined by the black hole spin.In the simulations with  * = 0.9, there are some indications of power around  max at later times, although nothing significant compared to the features observed at intermediate frequencies.On the other hand, in the early phase of the  * = 0.5 simulation, there is strong evidence of power along  max , though it is challenging to discern whether this power originates from trapped g-modes or is associated with vertical epicyclic motions.The thindisk dispersion relation (eq.14 in Ferreira & Ogilvie (2008) for the isothermal case) predicts  = 0 inertial wave propagation where the squared, Doppler-shifted wave frequency is less than the squared, radial epicyclic frequency, i/.e.,  2 <  2 .However, Fig. 8, particularly for simulation a5b15L3, shows that the power near  max occurs at frequencies greater than the local epicyclic frequency.Therefore, while the trapped g-modes are a possibility in tilted disc simulations, their presence is actively disfavoured by these PSDs.

QPO features
We now shift our focus to the prominent QPO-like features present at intermediate frequencies.The feature at 0.002  3 /  (shown in Fig. 8 & 9), which is present in all the tilted simulations, seems to correspond to the radial/vertical epicyclic frequency near 20  / 2 .This frequency is associated with the outer edge of the precessing torus, right around where the Lense-Thirring timescale becomes shorter than the accretion timescale (see Fig. 4).It is intriguing to see that the torus, while precessing as a rigid body within  ∼ 20  / 2 , also undergoes global epicyclic oscillations at a frequency corresponding to its outer edge.The cyan dashed and dotted horizontal lines around 0.002  3 /  in Fig. 9 represent the radial and vertical epicyclic frequencies at 19  / 2 , respectively.
The relatively high-frequency feature in these simulations is more complex.In the low-resolution simulations, the feature at around 0.005  3 /  matches with the vertical epicyclic frequency at the location of the initial pressure maximum of the torus, which is Figure 8. Power-spectral densities for all four simulations computed for the density-weighted average of the polar velocity,   .The solid, dashed, and dotted curves correspond to the Keplerian, radial, and vertical epicyclic frequencies, respectively.We note power along or between the epicyclic frequencies in all simulations, as well as horizontal (QPO-like) features in the tilted simulations.
Figure 9. Dynamical power spectra of ⟨  ⟩  ,  ( ,  ) at  = 2  / 2 for all three tilted disc simulations, computed using a sliding window of length 3000  / 3 .The cyan, dashed and dotted lines correspond to the radial and vertical epicyclic frequencies, respectively, at 19  / 2 , while the cyan dash-dot line represents the sum of the Keplerian and radial epicyclic frequencies.The long-dashed black lines represent the maximum radial epicyclic frequencies, and the black, dashed and dotted lines in the second and third panel correspond to the radial and vertical epicyclic frequencies, respectively, at the radius of the initial pressure maximum, i.e., 9  / 2 .9  / 2 .This frequency is shown by the black, dotted horizontal line in Fig. 9.However, this correspondence does not seem to hold for the high-resolution simulation.In a9b15L4, the higher frequency feature is instead observed around 0.004  3 / , which coincides with the sum of the Keplerian and radial epicyclic frequencies at 19  / 2 (indicated by the cyan dash-dot line).Our current understanding is that the mesh refinement at 10.4  / 2 may be influencing this feature, potentially enhancing it over the epicyclic frequencies at 9  / 2 , given that the local vertical epicyclic frequency at the mesh refinement radius happens to also be close to 0.004  3 / .
In the literature, tori examined for the study of various oscillation modes are typically isolated and exhibit polytropic behavior.The oscillation frequencies of such tori are closely tied to the characteristic frequencies at their pressure maxima.This connection arises because even if each annulus of the torus oscillates at its local epicyclic frequency, the average oscillation frequency converges toward the epicyclic frequency of the pressure maximum, given that most of the torus's mass is concentrated there.Consequently, it is unsurprising to observe that the feature at 0.002  3 /  corresponds to the radial/vertical epicyclic frequency at approximately 19  / 2 , as that is where the bulk of the mass resides in our simulations.Although we initialized the torus with the pressure and density maximum at 9  / 2 , the presence of the outer thin disc, with a density at least an order of magnitude higher, causes the pressure maximum and density maximum of the precessing torus to migrate towards its outer edge (see Fig. B1).Additionally, the strong coupling of the thick disc within 20  / 2 , induced by the bending waves due to the Lense-Thirring effect, further encourages global oscillations.It is thus intriguing to observe the profound influence of the outer thin disc on the inner torus, as it modifies its shape and mass distribution, and consequently, promotes eigenfrequencies that diverge significantly from those of an isolated torus.
Furthermore, we have observed that the Brunt-Väisälä frequency exceeds the local epicyclic frequencies within the torus region but is comparable to the local epicyclic frequencies in the outer thin disc region.Therefore, it is plausible that the oscillations in our simulations are influenced to some extent by a buoyancy restoring force.While this may not significantly alter the types of nonlinear interactions that could excite oscillations, it may complicate comparisons with oscillations in isolated tori, as seen in previous studies (e.g., Blaes et al. 2006;Mishra et al. 2017).

Radial and Vertical epicyclic motions
The presence of power over a range of radii at a fixed frequency suggests that the torus is undergoing the motions of some global mode.Earlier studies (e.g., Fragile & Blaes 2008) reported global radial epicyclic motions in tilted disc simulations, but those motions were not correlated with any particular QPO frequency.To further investigate the nature of these oscillations, we compute the PSD of the density-weighted average velocities for different modes, represented as ⟨  cos() sin( −  mid )⟩  ,  /⟨⟩  ,  , for  =  and .Here,  represents the azimuthal wavenumber and  denotes the number of nodes in the vertical direction7 .For this analysis, we compute the eigenfunctions in the transformed coordinates, where the midplane of the disc is aligned with the  ′ - ′ coordinate plane so that  ′ mid = /2.
The results remain the same if we instead repeat the analysis in the coordinate frame using  mid from Eq. 13.Fig. 10 displays the PSD of   ′ for all combinations of  and  between 0 and 1.The overlaid curves carry the same meaning as in Fig. 8.In addition to these curves, we include a dash-dot curve representing Ω K +  and a vertical solid line at 19  / 2 as a reference point.Here Ω K and  represent the Keplerian and radial epicyclic frequency, respectively.We find power in the  = 1,  = 0 panel that resembles radial epicyclic motion, e.g., the entire torus oscillating around the black hole radially.The power in  = 0,  = 1 maybe associated with the radial pulsations with a node at the disc midplane, where one half of the disc is radially expanding while the other half is radially contracting.Additionally, there is power along the Ω  +  curves in the outer, thin-disc region, indicating that while the thin disc oscillates at the local epicyclic frequency, the thick disc undergoes global oscillations.
Fig 11 shows the same analysis as in Fig. 10, but for   ′ .Note that we now observe power in the opposite panels, namely  = 0,  = 0 and  = 1,  = 1.The  = 0,  = 0 mode indicates pure vertical oscillations in the disc, while the  = 1,  = 1 mode corresponds to a motion resembling non-axisymmetric, vertically pulsating mode.The power spectra, including mode numbers > 1, reveal that combinations with odd differences between the mode numbers show power in   ′ , while combinations with even differences between  and  exhibit power in   ′ .The oscillations in   ′ and   ′ are also prominently evident in the respective space-time diagrams, as shown in the Appendix (see Fig. B2).
The excitation mechanism for these oscillations and the nature of these modes are not entirely clear.However, previous analytical calculations (Kato 2004(Kato , 2008;;Ferreira & Ogilvie 2008) and subsequent simulation studies (Henisey et al. 2009;Dewberry et al. 2020b) have shown that disc warps and eccentricities can excite and amplify inertial modes in the disc through non-linear coupling.These excitation mechanisms can be described as nonlinear three-wave couplings between a nearly zero-frequency distortion with  ∼ 1, and a pair of waves with nearly equal frequencies and  values differing by one.The appearance of both  = 0 and  = 1 power at the same frequencies in Figs. 10 and 11 indicates that such a parametric-type instability may be responsible for the oscillations in our simulations.

Vertical oscillations
Even without the aid of a PSD, the vertical oscillations reported above can be clearly seen in a plot that tracks the disc midplane angle,  mid , over time, as we show in the top panel of Fig 12 .We compute  mid based on the concentration of the disc's mass: For clarity, we show the variations of  mid over time for a specific azimuthal angle ( = 0) at different radii.At  = 5  / 2 , clear oscillations in  mid can be observed in each case, indicating an upand-down motion of the disc.The amplitude of these oscillations decreases as one moves away from the black hole.The bottom panel of Fig. 12 provides insight into the azimuthal distribution and evolution of these vertical oscillations over time.A dominant  = 1 pattern is evident, arising from the deviation of the disc midplane from the grid plane  = 0.Although the disc starts with  mid = /2 for all , once the precession begins, half of the disc has  mid < /2 while the other half has  mid > /2.The  = 1 pattern shifts with time due to the precession of the thick disc.The black, dashed curve in Fig. 12 corresponds to the precession rate evaluated at 5  / 2 obtained directly from the simulation data.In addition to the evidence this plot provides for tilt and precession, we also observe additional regular and moderately rapid, oscillations of the maximum and minimum values of  mid , corresponding to the oscillations seen in the top panel.The period of these oscillations is approximately 500  / 3 , matching the frequency of 0.002  3 /  noted in the PSD of polar velocity (Fig. 8).In other words, we are seeing the global  = 0,  = 0 vertical oscillations discussed earlier.

Comparison with previous works
In this paper, we have investigated the timing and dynamical properties of tilted, truncated accretion discs.While the truncated disc picture is rather strongly motivated by observations of the hard state in X-ray binary systems, their formation from first principles is not yet fully understood.Recent simulations have made progress in selfconsistently forming truncated discs by incorporating strong magnetic fluxes.In these simulations, a thin disc initiated with strong magnetic flux results in the accreting gas dragging in enough poloidal magnetic flux to saturate the black hole and disrupt the accretion flow.This leads to the formation of a thin disc that is truncated beyond the normal marginally stable radius, with a hot, magnetically saturated corona present inside (Mishra et al. 2022;Liska et al. 2022;Dihingia et al. 2022).Another proposed mechanism for truncated disc formation, known as the 'disc evaporation model' (Meyer & Meyer-Hofmeister 1994;Liu et al. 1999), suggests that the inner regions of the disc are evaporated at low densities.However, to further investigate phenomena like disc evaporation, additional physics such as thermal conduction needs to be incorporated into simulations, which would be both challenging and computationally expensive.
Numerical simulations conducted thus far suggest that strongly magnetized discs do not undergo precession (McKinney et al. 2013;Ressler et al. 2023;Fragile et al. 2023).This is attributed to the magnetic torque exerted by the strong fields accumulated near the black hole, which effectively aligns the disc with the symmetry plane of the black hole.Since our primary motive is to study precession in a truncated disc, we chose a magnetic field configuration that leads to a weakly magnetized corona in the inner regions (like a Standard and Normal Evolution (SANE) accretion flow Narayan et al. 2012) and maintain the outer, truncated thin disc with the aid of artificial cooling.In the future, we wish to explore the precession of strongly magnetized discs with a truncated disc geometry.
Geometrically thin discs with extreme tilt angles in the presence of a spinning black hole can experience disc tearing, as highlighted in previous studies (Nixon & King 2012;Raj et al. 2021;Liska et al. 2021).Notably, Musoke et al. ( 2023) conducted further analysis of the tearing disc simulations presented in Liska et al. (2021) and identified high-frequency quasi-periodic oscillation (HFQPO) features corresponding to the radial epicyclic frequencies at the tearing radius.However, it is crucial to emphasize that our simulations differ significantly from the disc tearing scenario.In our simulations, there is no occurrence of disc tearing, and the torus exhibits global oscillations across a range of radii.
Oscillations at epicyclic frequencies have previously been reported in previous global relativistic hydrodynamic simulations of pressuresupported (non-accreting) tori around black holes (Schnittman & Rezzolla 2006;Mishra et al. 2017).Yet, in GRMHD simulations of accreting tori, there has been a lack of compelling evidence supporting the existence of QPOs associated with these oscillations (Schnittman et al. 2006).Conversely, GRMHD simulations of tilted accreting tori initially reported global variability at specific orbital frequencies (Henisey et al. 2009), only to later identify this as a transient feature (Henisey et al. 2012).It was also suggested that the presence of shocks in tilted discs enhances this variability, thus explaining the limited variability seen in aligned, accreting tori.Our results agree with the latter aspect of these studies, as we do not observe any variability in the aligned case.However, in contrast to the earlier studies on tilted accreting tori (Henisey et al. 2009(Henisey et al. , 2012)), we have identified a persistent QPO feature at 0.002,  3 /  in all our tilted disc simulations8 .This feature appears to correspond to epicyclic frequencies at the outer edge of the precessing torus, influenced by the presence of the outer thin disc -a component absent in any previous simulations documented in the literature.
Our findings contradict the conclusions of Marcel & Neilsen (2021), who argued that hot, thick accretion flows must be supersonic, and therefore cannot precess as solid bodies.Our results, as depicted in Fig. 4, demonstrate that our accretion discs are not supersonic, even in the hot, thick regions.In fact,   <   except for the very innermost regions.Furthermore, both in our previous work (Bollimpalli et al. 2023) and here (see Fig. 7), we demonstrated that the hot, thick disc precesses as a rigid body, albeit at a rate slower than predicted for an isolated disc, owing to the presence of the outer thin disc.The only question that remains is whether our simulations can reproduce the spectral and timing properties of the luminous, hard state, particularly the prominent type-C QPO.We plan to address this question in future investigations.

Observational relevance
Often, models attempting to explain QPOs as being associated with the epicyclic frequencies in geometrically thin discs face challenges in matching the observed frequencies due to the radial dependence of these frequencies.For example, HFQPOs that are commonly found in 3:2 frequency ratio pairs (Remillard & McClintock 2006) are suggested to be due to a resonance (Kluzniak & Abramowicz 2001a;Abramowicz & Kluźniak 2001), possibly a parametric resonance (Kluzniak & Abramowicz 2002) of oscillations at the radial and vertical epicyclic frequencies.In contrast, global oscillations in geometrically thick discs are also a viable explanation, as suggested by various studies (Rezzolla et al. 2003;Schnittman & Rezzolla 2006;Blaes et al. 2006;Mishra et al. 2017).Our findings, as discussed in Sections 4.2 and 4.3, provide compelling evidence that global oscillations can only be excited in the thick disc occupying the inner region of truncated discs, and only when the discs are misaligned.Thus, our results strengthen the argument that the presence of tilt is a crucial requirement for the occurrence of QPOs.
In addition to the low-frequency QPO, the PSD of black hole transients reveals a more complex structure, typically analyzed using a combination of Lorentzian components (Belloni et al. 2002).These components yield characteristic frequencies, but due to their broad nature, determining a single characteristic frequency, as is the case with a narrow peak, becomes a model-dependent task.During the early phases of the state transition, it is observed that as the source flux increases, these characteristic frequencies increase (Belloni & Stella 2014;Bhargava et al. 2021).This increase is often attributed to the inward movement of the truncation radius.Our findings, which associate the QPO feature at  ≈ 0.002  3 /  with epicyclic frequencies at the transition radius (i.e., the outer edge of the precessing torus), are in line with this understanding, since the epicyclic frequency increases with decreasing radius.To further validate this connection, we intend to conduct future truncated disc simulations with different torus sizes, representing varying transition radii, to examine their impact on QPO frequencies.
As per the findings of Fragile et al. (2016), the correlation between the QPO frequencies observed in low mass X-ray binaries could be explained by associating the type-C QPO with the torus precession (equivalent to the  = −1 vertical epicyclic mode), and the HFQPOs, often occurring in a near frequency 3:2 ratio, with the axisymmetric modes of the torus, specifically the breathing and vertical epicyclic modes, respectively.The current results from our simulations support this model, where we clearly see a  = 0 vertically epicyclic mode for the lower high frequency QPO in addition to the precession of the torus, which could be interpreted as an  = −1 vertical epicyclic mode.Since our simulations do not cover a full precession cycle, it is not possible to observe any such mode in our PSDs.We note that the other HFQPO we found in our simulations is not a pure breathing mode as in the Fragile et al. (2016) model, but rather something similar to a non-axisymmetric, vertically pulsating mode.
Another point is that HFQPOs are only observed in highluminosity states at specific hardness ratios (McClintock & Remillard 2006;Belloni et al. 2012).While this could be a selection effect due to the low signal-to-noise ratio in the low-luminosity states, it could also be that QPOs produced from global oscillations only have a detectable quality factor for small sizes of tori, reached as the source transits to the soft state.
Further, recent multiwavelength observations of black hole X-ray binaries have enabled detection of LFQPOs both in IR/Optical and very high energies(> 200 keV), at nearly the same frequency as the type-C QPO observed in X-rays (e.g.Kalamkar et al. 2016;Ma et al. 2021).The popular explanation for such multi-wavelength observation of a common QPO frequency is a simultaneous precession of the jet along with the corona (Malzac et al. 2018).Additional correlated QPO frequencies may be observed in the future if, for instance, the vertical oscillations observed in the hot, thick inner disc of our simulations were to propagate into the jet, which may happen if this region of the disc plays a role in confining or collimating the base of the jet.

SUMMARY
In this paper, we performed a suite of GRMHD simulations of truncated accretion discs, tilted and untilted, to study the dynamics and precession of the hot, inner, geometrically thick flow in the presence of a cold, outer, geometrically thin disc.Our main conclusions are: • There is a slow down of the precession rate: The precession rate of the inner torus in a truncated geometry is significantly slower than the precession rate of an isolated torus of the same size, as we noted previously in Bollimpalli et al. (2023).This decrease in the precession rate is caused by the flux of angular momentum between the two disc components.
• Tilted discs exhibit significantly more variability: We observe QPO-like features at different frequencies in all of our tilted disc simulations, whereas such features are absent in the untilted case.
• The inner, misaligned torus exhibits coherent vertical oscillations and radial epicyclic motions: We identify a prominent QPO-like feature at a frequency of approximately 0.002  3 / , which is consistently present in the PSDs of all three tilted simulations.This frequency corresponds to approximately 40 Hz for a black hole of 10  ⊙ .This feature seems to correspond to the radial or vertical epicyclic frequency at the transition radius.
We also observe another QPO-like feature at approximately 0.004  3 /  in the high-resolution simulation, a9b15L4, and approximately 0.005  3 /  in the low-resolution simulations, a9b15L3 and a5b15L3; in the a5b15L3 case its frequency could be in a 3:1 ratio to the fundamental.We hypothesize that the feature in the high-resolution simulation corresponds to the vertical epicyclic frequency at the mesh refinement radius, while the feature in the low-resolution simulations corresponds to the vertical epicyclic frequency at the radius of the initial pressure maximum in the torus.
• The thin disc exhibits variability at the local epicyclic frequency corresponding to each radius: In addition to the QPO-like features associated with the radial and vertical epicyclic motions of the inner thick disc, our simulations reveal power along the Keplerian frequency curve in the thin disc region, indicating local oscillations at the corresponding Keplerian frequencies.
Overall, our results align well with the observations of HFQPOs in the black hole X-ray binary systems.

APPENDIX A: SIMULATION OF AN ISOLATED TORUS
To ensure that the slowdown in the precession rate is a physical effect caused by the outer thin disc rather than a numerical artifact, we performed one additional simulation of an isolated torus using the exact same numerical setup as in simulation a9b15L4, but excluding the thin disc beyond 15  / 2 .Fig A1 shows the time evolution of the precession angle.We find that the torus completes one full precession cycle within 4500  / 3 .This is consistent with the expected precession rate for an isolated torus extending from 5 to 15  / 2 (see eq. 2 in Ingram et al. 2009), but much faster than what we see when the thin disc is included.

APPENDIX B: SUPPLEMENTARY FIGURES
Fig. B1 shows radial profiles of the time-averaged surface density (in thicker curves) compared to their initial profiles (in thinner curves) for all three tilted disc simulations.The time averaging is done from 3000  / 3 to the end of each simulation.Note that the torus in the initial profile has a density maximum at 9  / 2 independent of the thin disc whose surface density is at least an order of magnitude higher than the maximum torus density.Due to this, the maximum density of the torus migrates towards the outer thin disc as the simulation evolves.
To further substantiate the presence of the reported radial and vertical oscillations in Section 4.3, we present spacetime diagrams of the azimuthally and vertically averaged   ′ (top panels) and   ′ (bottom panels) for simulations a9b15L4 and a5b15L3 in Fig. B2.To emphasize the  = 0,  = 0 modes observed in Fig. 11, we employ azimuthal averaging over the entire 2 and vertical averaging over the disc spanning ±25 • around the midplane.For   ′ , vertical averaging considers only the upper half of the disc, specifically  ∈ [65 • , 90 • ].The oscillations around ∼ 0.002,  3 /  are clearly discernible in simulation a5b15L3, while oscillations at relatively higher frequencies (∼ 0.004,  3 / ) are more pronounced in simulation a9b15L4, particularly in   ′ , at least for the time period depicted in this plot.

Figure 2 .
Figure 2. Space-time plot of the disc scale height, /, for the 4-level untilted (a9b0L4, top) and tilted (a9b15L4, bottom) disc simulations.The two-flow geometry is maintained well throughout the evolution in both cases.

Figure 3 .
Figure 3. Top: Space-time plot of the mass accretion rate for simulation a9b15L4, measured in arbitrary units.Bottom: Evolution of the mass accretion rate at two radii, 5 and 25  / 2 , representing the thick-and thin-disc regions, for all four simulations.Note the differing scales of the two plots.

Figure 4 .
Figure 4. Characteristic timescales (time-averaged over the interval  ∈ [10000  / 3 ,  end ]), scaled to the local orbital period, for all three tilted disc simulations.The Lense-Thirring (LT) precession timescales are shorter than the accretion timescales in the thick disc region (6 ≲ / g ≲ 20), while the sound-crossing timescales are shorter than the precession timescales everywhere except in the very inner regions of the thick discs, thus allowing for rigid-body precession of the thick disc.

Figure 6 .
Figure 6.Space-time plot for the tilt angle,  (measured in degrees), for the a9b15L4 simulation.The blue-solid curve in the top panel shows the radial profile time averaged over 10000 − 25000  / 3 .

Figure 7 .
Figure 7. Space-time plot for the rate of change of the precession angle.The bending wave, which propagates at roughly half the sound speed, is depicted by the black, dashed curve.The horizontal and the vertical dotted lines mark the time 1300  / 3 and the truncation radius 15  / 2 , respectively.The consistent colour in the precession rate within the range 6 − 15  / 2 suggests rigid-body precession of the thick-disc region.

Figure 12 .
Figure12.Top: Time evolution of the disc's midplane angle,  mid , measured at a fixed azimuthal angle,  = 0, for radii  = 5, 10, 15, and 20  / 2 for the three tilted-disc simulations.Bottom: A time-space plot of the midplane angle,  mid , at  = 5  / 2 for the simulation a9b15L4.A clear  = 1 pattern is seen which precesses along with the rest of the thick-disc region (black, dashed line) and also oscillates with a period close to 500 g /.

Figure A1 .
Figure A1.Space-time evolution of the precession angle,  (measured in degrees), for an isolated torus simulation.

Figure B2 .
Figure B2.Space-time diagrams showing the azimuthally and vertically averaged radial (  ′ , in top panels) and poloidal (  ′ , in bottom panels) velocities for simulations a9b15L4 (left column) and a5b15L3 (right column).To illustrate the  = 0,  = 0 modes in   ′ , we consider vertical averaging over ±25 • around the mid-plane, while for the  = 0,  = 1 modes in   ′ , we consider only the upper half of the disc.

Table 1 .
Simulation parameters luminous hard state