X-ray Quasi-periodic Eruptions driven by Star-Disc Collisions : Application to GSN069 and Probing the Spin of Massive Black Holes

X-ray quasi-periodic eruptions (QPEs) are discovered recently in active galaxies with unknown driven mechanism. Under the assumption that QPEs are caused by star-disc collisions, we adopt full relativistic method and show that both the orbital parameters of the star and also the mass and spinning of the massive black hole (MBH) can be revealed by using the time of arrival (TOA) of the QPEs. By applying the model to the observed QPEs of GSN069, we find that the star is in a near-circular orbit ( $e_\bullet=0.05^{+0.02}_{-0.02}$) with semimajor axis of $\sim 365^{+54}_{-49}r_{\rm g}$ around a MBH with $M_\bullet=3.0^{+0.9}_{-0.6} \times10^5M_\odot$. The alternative short and long recurring time of the QPEs of GSN069 can be well explained by the small eccentricity and the orbital precession of the star. We find that the QPEs of GSN069 are possibly driven by a striped stellar core colliding with accretion disc after partial tidal disruption event around the MBH. For GSN069-like galaxies, if continuous X-ray monitoring of QPE events can be accumulated with uncertainties of TOA $\lesssim 100-150$s, the spin of massive black hole can be constrained by fitting to QPEs. Our results show that the timing of QPEs can provide a unique probe for measuring the spinning of MBH and tests of no-hair theorem.


INTRODUCTION
X-ray quasi-periodic eruptions (QPE) are rapid and extremely intensive bursts of X-ray emission which repeat about every few hours from regions around massive black holes (MBH) in galactic nuclei. QPEs have been found recently from the active galactic nucleus of GSN069 (Miniutti et al. 2019) and galaxy J1301.9+2747 (Giustini et al. 2020). More recently, QPEs are detected in the nucleus of two previously quiescent galaxies (Arcodia et al. 2021). Multiple epoch of observations suggest that QPEs may last from at least months to decades (Miniutti et al. 2019;Giustini et al. 2020). Currently it is still unclear of the underlying physical mechanism that drive the burst of these events. Some of the possibilities are the limit-cycle oscillations induced by instabilities of the accretion flow (Lightman & Eardley 1974) or an edge-on binary black hole gravitationally lensing the light from each others' accretion disc (AD) (Ingram et al. 2021). However, they are both found inconsistent with the properties of the eruptions (Arcodia et al. 2021).
Other possiblities are the falling clumps of partial TDE (Coughlin & Nixon 2020), or the outburst arised due to the misalignment between the disc and the spin of the MBH (Raj & Nixon 2021). One of the tempting mechanisms is that the Xray QPEs are driven by an orbitor around the MBH that is with masses much smaller than the MBH (Arcodia et al. 2021). One of such examples is a white dwarf partially disrupted by MBH when crossing the periapse per orbit (King 2020).
It has been suggested that quasi-periodic flares may be produced by collisions between a star and an accretion disc surrounding MBH (Karas & Vokrouhlicky 1994;Nayakshin et al. 2004;Dai et al. 2010). Due to collision on the disc, gas are shocked and subsequently a bright hot spot may appear and evolve on the disc (Zentsova 1983;Ivanov et al. 1998;Nayakshin et al. 2004;Pihajoki 2016;Semerák et al. 1999;Pasham et al. 2019). The star-disc collision model may be able to explain the X-ray periodicity of galaxy RE J1034+396 (Dai et al. 2010;Gierliński et al. 2008), and the optical variability of source OJ287 (Lehto & Valtonen 1996;Dai et al. 2010). By adopting post-Newtonian methods, Valtonen et al. (2008Valtonen et al. ( , 2011 show that the timing of these flares may be a useful probe for general relativity in strong field. Here we study the timing of X-ray QPEs under the assumption that they are driven by the star-disc collisions. We use a full general relativistic numerical method (developed in Zhang et al. 2015;Zhang & Saha 2017) to simulate the orbit of a star and also the propagation of lights emitted from the disc to the observer. We then develop a numerical method which can extract the orbital parameters of the star from the timing of X-ray QPEs. We apply the method to GSN069 and set constraints the orbit of the star. We also discuss the possiblity that the orbitor in GSN069 is a result of partial tidal disruption event. We further show that using QPE alone can set strong constraints on the spinning of the massive black hole if QPEs can be observed with high enough timing accuracy. Details of the modeling and the results are shown in the following sections.

MODEL AND METHODOLOGY
The prime assumption is that the X-ray eruptions are due to super-sonic collisions between a star 1 bound to the MBH and a standard geometrically thin AD, similar to those in Rauch (1995) and Dai et al. (2010). Here we adopt full general relativistic numerical framework developed in Zhang et al. (2015) (See also Zhang & Saha (2017)) to calculate the orbit of the star and the propagation of lights from the disc to observer under Boyer-Lindquist coordinates.
Suppose that the peak of each QPE corresponds to the time when the star intersects with the midplane of the disc.
At each time of intersection, we use raytracing technique to trace back a number of photons from a distant observer until we find the one that hits on the position within distance of < 10 −8 r • from the star, where r • is the distance of the star to the MBH. The ray-tracing is fast and accurate as we use Jacobian elliptic functions and Gauss-Kronrod integration scheme 1 We will show later that the orbitor is actually more consistent with a remnant core of a red giant after partial tidal disruption event. Nevertheless, we generally call the orbitor "a star" across the paper.
for the integration of motion equation of photons (For more details see Zhang et al. (2015)). Then the time of arrival (TOA) t TOA of the eruption due to collision measured in the observer's frame is then given by where t is the coordinate time of the star when it intersects with the disc mid-plane and t prop is the time of propagation of the photon from the intersection to the observer found by the ray-tracing method. For simplifity, we ignore the secondary image of the eruption produced by photons running around the other side of the MBH and twisting back to the observer, as usually the amplitude of these high order images are small (Karas & Vokrouhlicky 1994;Dai et al. 2010). We also ignore the gravitational wave orbital decay and the possible deviations of the star orbit due to collision for now, as we find later that they can be safely ignored for GSN069 (See more details in Section 4). Figure 1 illustrate examples of trajectories of orbitors and photons in Newtonian case (left panel) and relativistic cases under Schwarzschild and spinning MBHs (middle and right panel). We can see that the orbit of the star intersects with the accretion disc twice per orbit, with the time of intersection separating with alternative short and long intervals if there is a non-zero orbital eccentricity. In the meanwhile, the orbit of the star is precessing due to relativistic effects (Schwarzschild orbital precession and etc), thus the location of the intersection will change in each revolution. Due to the curved spacetime under Schwarzschild or Kerr metric, the propagating time of photons from the position of the intersection to the observer is also different. Combining these effects, the TOA of each eruption may then appear with irregular periods.
As the TOA of X-ray QPEs encode various orbital and relativistic effects, in principle the parameters of the MBH or orbit of the star should be reconstructed by the timing of QPEs. We then apply our method to the QPEs observed in GSN069, with more details shown in the following section. As we find that the current observations of GSN069 are not able to recover the spin parameters of the MBH, we discuss the constraints of spin parameters from timing of QPEs for future observations in Section 5.

EXTRACTING THE ORBIT OF STAR BY X-RAY TIMING OF QPES IN GSN069
Multiple QPEs on Seyfert galaxy GSN 069 have been reported by (Miniutti et al. 2019), which recurrent every ∼ 9hr with duration of about ∼ 1hr around a MBH of . For illustration purpose, the observer is located at distance of r = 30rg from the MBH (for simulations performed in this work, it is actually at r = 10 7 rg). Note that the ascending node or descending node of the orbitor (the red solid line) remain constant for a Schwarzschild MBH, but change in each revolution for a spinning MBH due to the precession of the orbital plane of the orbitor.

166.5
Note-a The prior ranges for parameters. b The constraints of parameters from mean-likelihood in 2σ confidence level. c The parameters with minium χ 2 value (χ 2 = 16.2). d Here α• should avoid angle near 0 or 180 degree where the disc is edge on, which can cause difficulties in the light tracing numerical methods. In principle, the MBH-star system contains in total of 9 independent parameters: 3 independent parameters for MBH: mass (M • ), the dimensionless magnitude of spin (a), the line of sight inclination of spin (i, or the normal direction of accretion disc if a = 0); 6 parameters of the orbit of star that defined on the disc plane: the semimajor axis (SMA, a • ), eccentricity (e • ), inclination (I • ), ascending node (Ω • , defined respect to the projection line of sight on the disc plane.), argument of periapsis ω • , and the true anomaly f • . Parameter M • can also be replaced by the Keplerian orbital period Due to the limited timing accuracy of each TOA, we found that the spin of MBH can not be reconstructed for GSN069, thus we fix a = 0. If alternatively set a = 0.99, we find that the constraints on other parameters are only weakly affected. The independent parameters of the model in case of a Schwarzschild black hole reduce to 6 (Dai et al. 2010 where α • is the angle between the vector pointing from the MBH to the ascending node of the star and the vector pointing from the MBH to the observer, and cos α • = sin i cos Ω • . We then use the Markcov-Chain Monte-Carlo (MCMC) method to reconstruct these parameters. As the observations are on discrete time periods with different number of eruptions, it is difficult to correspond an eruption in observation to that in simulation. Thus, here we evaluate the χ 2 value by a method similar to those in the pulsar timing. Suppose that there are a sequence of QPEs indexed by k = 1, 2, · · · and the TOA of the k-th one is given by t TOA,k = F(k). In observation, not all of them are covered. If there are N discrete observations, and the number of QPEs in each one is given by M i , i = 1, · · · , N , then χ 2 is given by where    index k. I j=1 is an integer which is closest to the index given by F −1 (t TOA,i1 ). σ tij is the measurement error of the TOA of the j-th flare of the i-th observation. Due to the random brightness fluctuations of the QPE events, it is difficult to get directly the peak of each QPE from the observed lightcurves. We adopt the Bayesian block method (Scargle et al. 2013;Vaughan 2010) which can split the time series to optimal segments to help to identify the peak location of each QPE event. The analysis results for the observed QPEs of GSN069 are shown in Figure 2. The TOA of each QPE is reasonably expected to be located at the center of the Bayesian blocks with maximum flux of a local time series, and the 1σ measurement error of TOA is assumed the half of the time span of that block. The uncertainty of the TOA for QPEs in GSN069 are found between 430−875s.
The reconstructed parameters are shown in Figure 3 and Table 1. We find that the TOA of QPEs of GSN069 are consistent with an orbiter with SMA of ∼ 365r g and a small eccentricity of ∼ 0.05 (according to the mean likelihood, marked by light-blue cross in each panel of Figure 3). The angle α • ∼ 106 • suggest that the light of sight direction is tilted to nearly perpendicular to the vector pointing from the MBH to the ascending node of the star intersect with the accretion disc.
The best-fit value deviating slightly from those from mean likelihood in Table 1, is due to two local minimum of χ 2 , as more clearly shown in Figure 3. The two local minimum is more apparent for M • and a • , where one of them are the best-fit value at M • = 2.69 × 10 5 M and a • = 390r g (marked by yellow cross), and the other are at M • ∼ 3.8 × 10 5 M and a • ∼ 320r g . As the size of the MCMC chains is sufficiently large (∼ 10 5 accepted MC samples), the two local minimum is more likely due to the fluctuations caused by small sample size of the QPE or some unknown systematics in the data.
Comparison of the observed QPEs of GSN069 and the TOA of each eruption in the best-fit model are shown in Figure 2 (marked by red arrows). The TOA of each QPE from observations are well consistent with those from simulations within 2σ errors.
The time intervals between QPEs in each observation are more clearly presented in Figure 4. It is interesting to see that the TOA of both the observed and simulated QPEs are with alternative short and long intervals, which irregularly vary by amount of 2 − 3ks. Such quasi-periodic behaviour can be explained by combination effects of the Schwarzschild precession and the presence of a small eccentricity of the orbit of the star. The time interval between QPEs should be in orders of P 1/2 + δP (f 0 ) where P 1/2 is the half orbital period and δP (f 0 ) ∼ 4e sin f 0 /ω • = 4e sin f 0 (M • G) −1/2 a 3/2 • (when e • 1), whereω • = (GM • /a 3 • ) 1/2 is the orbital angular velocity, f 0 is the true anomaly of the orbitor when it hits on the disc. If there is no orbital precession (f 0 = constant), the time interval should be either P 1/2 +δP (f 0 ) or P 1/2 +δP (f 0 +π) = P 1/2 −δP (f 0 ). Due to Schwarzschild and Lense-thirring orbital precession, f 0 precess in each orbit and the interval of QPEs will be modulated between P 1/2 + δP (π/2) and P 1/2 − δP (π/2). For GSN069, we find that δP (π/2) ∼ 2.4ks, which is consistent with those shown in the right panel of Figure 4. In Schwrzschild black hole, the orbital precession become 2π after in total of N p,S flares, which is where δω S = 6π(a • /r g ) −1 (1 − e 2 • ) −1 . This is also consistent with the period of modulation N p,S ∼ 258 shown in the right panel of Figure 4. The difference of time per flare is then approximately 2δP (π/2)/N p,S ∝ a 1/2 • 2 , thus the larger the distance of the orbitor, the larger the relativistic effects on TOA of QPEs around a Schwarzchild MBH.
As a summary, we find that the observed TOA of the QPEs of GSN069 can be well explained by the star-disc collisions models and some of the parameters of the orbit of the star can be well constrained by using only the timing of the X-ray emission of QPEs. The quasiperiodic behaviour can be explained by combination effects of a small orbital eccentricity of the orbitor and the Schwarzschild precession.

IS THE ORBITOR IN GSN069 A STELLAR CORE FROM A PREVIOUS PARTIAL TIDAL DISRUPTION EVENT?
Recent observations suggest that GSN069 exhibit some characteristic of partial tidal disruption event (TDE) (Zhao et al. 2021;Sheng et al. 2021). Considering that the QPEs of GSN069 exhibit long-lived signatures of tidal disruption event (Shu et al. 2018), here we show that, if the star-disc collision model is the driven mechanism of QPEs of GSN069, the orbiter may be a stellar core remained after a partial tidal disruption of an evolved star (possibly a red giant) around the MBH in GSN069. More details are provided as follows.
Partial TDEs can happen if a red giant star pass through the MBH near the tidal radius r p ∼ r t = r (M • /m ) 1/3 , where r and m is the radius and the mass of the red giant (e.g., Bogdanović et al. 2014;Coughlin & Nixon 2019;Chen & Shen 2021;MacLeod et al. 2012). Typically for red giant r = 5 ∼ 20R when m ∼ 1 − 5M (MacLeod et al. 2012). As the stellar core is more compact than the outer envelope, the stellar core remained while the envelope is disrupted by the tidal force of the MBH. Then the stellar core should be initially with orbital periapsis r p ∼ 300 − 2075r g for GSN069 with M • ∼ 3 × 10 5 M . This result is consistent with the distance of ∼ 400r g of the orbitor from the MBH for GSN069 resulting from our constraints. Note that the above estimation is for initial value of r p and the orbit of the stellar core may evolved later due to the collisions with the accretion disc (See text later).
During each collision, part of the kinematic energy of the star are transferred into the shocked gas. As the radius of the remnant core is usually r core ∼ r /20 ∼ 0.1R (MacLeod et al. 2012), the X-ray luminosity of the eruption is then given by where f X is the fraction of emission in X-ray that depends on the details of radiative processes (Nayakshin et al. 2004), r core , α, and λ is the radius of the stellar core, viscosity, radiative efficiency and Eddington ratio of the disc, respectively. The typical X-ray luminosity of QPEs are in orders of 10 41 ∼ 10 42 erg/s (Miniutti et al. 2019;Arcodia et al. 2021), thus the above estimation is consistent with the observed luminosity when f X ∼ 0.01 − 0.1.
The star will experience drag force caused by sweeping up of the intercepted disc gas in each collision, and finally the orbit is circularized and aligned with the disc. Following Rauch (1995), the alignment timescale (t align ) is assumed approximately the time needed for star to intercept with the disc material with masses equal to its own masses, which is given by 3 t align t orb ∼ 1.3 × 10 6 m core 0.1M where t orb is the orbital period, m core ∼ 0.3 − 0.6M (MacLeod et al. 2012) is the mass of the stellar core to the MBH, respectively. For GSN069 with M • ∼ 3 × 10 5 M and a • ∼ 400r g (as shown in Section 3), it is then suggest that the QPEs can last for about ∼ 1170 yrs. The corresponding gravitational wave orbital decay is with timescale of ≥ 2 × 10 6 yr if e ≤ 0.6 and m core ∼ 0.3M . As the observation time of QPEs are typically in about 1yr, it is justify that in our simulations we ignore the changes of the orbit due to both the star-disc collision and the gravitational wave orbital decay.
As QPE sources are assumed results of partial TDE events, the event rate of QPE sources (R QPE ) is expected to be the same as those of partial TDE (R pTDE ), the later of which are expected to be about 10% of those normal TDE (disrupting a main-sequence star), i.e., R pTDE ∼ 0.1R TDE ∼ 10 −5 yr −1 per galaxy (with M • < 10 7 M ) (MacLeod et al. 2012). As the duration of QPE events can be T QPE ∼ 1000yr, the duty cycle of QPE sources is then expected to be around R QPE T QPE ∼ 10 −2 , suggesting that QPEs may appear in every one of a hundred late type galaxies (with M • < 10 7 M ) at any given moment. It will be more easy to understand the frequency of QPE sources if we obtain the number ratio between the QPE sources and those of normal TDE: where T TDE ∼ 0.1 − 1yr is the duration of normal TDE events. The above estimation suggest that theoretically the QPEs sources should be much more abundant than those of normal TDE sources. Currently there are only a few QPE sources revealed, possibly because the observation of them is quite expensive, and that not all of them produce bright enough flares.

CONSTRAINING THE SPIN PARAMETERS OF THE MBH BY TOA OF QPES
The TOA of QPEs can be also affected by the presence of the spinning of the MBH due to Lense-Thirring effects on the orbit and the additional deflection of lights by the frame-dragging effects. We investigate the constraints on the spin parameters of the MBH (also other model parameters simultaneously) by using mock observations, corresponds to which can possibly be collected in near-future for QPE sources. We find that the inclination I • has a near-degeneracy with the spin parameter a, which largely slows down the convergence of the MCMC simulations. As usually I • can only be poorly constrained, here we fix I • for simplicity. If I • is taken as a free parameter, we find that our following constraints on a is only slightly weaker (e.g., increase from ∆a ∼ 0.5 to ∆a ∼ 0.6). For demonstration purpose, we choose QPEs from a GSN069-like galaxy given the best-fit values suggested by Table 1: M • = 2.65 × 10 5 M (or P k = 63.16ks), a = 0.9, i = 60 • , a • = 390r g , e • = 0.06 or 0.5, There are in total of 8 free parameters in the MCMC: (P k , a, i, a • , e • , Ω • , ω • , f • ). We simulate in total of 425 flares, which corresponds to a time duration of ∼ 155 days. To mimic the intermittent observation, only the 1 − 5th, 106 − 111th, · · · and 421 − 425th eruptions (in total of 25 eruptions) are used for recovering.
By performing a number of MCMC simulations, we explore the constraints on the spin magnitude in case of different timing measurement uncertainties or the orbital eccentricities. The results are shown in Figure 5. In the case when e • = 0.06, the constraints of spin are only possible if σ t 100s. However, if e • = 0.5, the constraints can be stronger and the spin parameter can be revealed if σ t ≤ 150s. The results of the reconstructed model parameters are also shown in the left panels in Figure 5 given the uncertainties of TOA σ t = 50s and e • = 0.06. We can see that both the orbital elements and the spin of MBH can be recovered, where the uncertainties of ∆a ∼ 0.49 and ∆i ∼ 85 • , respectively.
In principle, orbiters similar to the above one but with higher eccentricities, or QPE observations covering longer time span, can set even stronger constraints on the spinning of MBH. By methods similar in Section 3, we expect that the time difference per flare due to the lense-thirring effect is approximately 2δP (π/2)/N p,L , where N p,L = 4π/|δω L |, and |δω L | = 12πa(a • /r g ) 3/2 (1− e 2 • ) −3/2 cos I • . Note that 2δP (π/2)/N p,L is independent of a • . Similarly, the variations due to quadruple momentum effects 2δP (π/2)/N p,Q ∝ a −1/2 • , where N p,Q ∝ a −5/2 • . We defer exploring more details of the constraining of the spin on various conditions of the orbitors by timing of QPEs to future studies.
6. SUMMARIES X-ray quasi-periodic eruptions (QPEs) in active galactic nuclei are found recently with unknown driven mechanisms. By assuming that QPEs are driven by star-disc collisions, we develop a full relativistic model to reconstruct the parameters of the orbiter and also the spinning of massive black hole (MBH) by using the time of arrival (TOA) of these QPEs.
We apply our method to the QPEs of GSN069, and find that the orbital parameters of the star can be well constrained, although the spinning of MBH can not, due to the limited timing accuracy of the observed samples. The orbiter in GSN069 is found with semimajor axis of about 365 +54 −49 r g and small eccentricity (e • = 0.05 +0.02 −0.02 ) around a MBH with M • = 3.0 +0.9 −0.6 × 10 5 M , We find that the variations on the recurring time of QPEs of GSN069, with repeating irregular alternative long and short intervals, can be well explained by the presence of the orbital eccentricity and the relativistic orbital precession. We also discuss the possibilities that the QPEs in GSN069 are caused by continual collisions between a stripped stellar core with the accretion disc around the MBH, the former of which are formed by partial tidal disruption of a red giant.
We find that, for a GSN069-like QPE source that last more than several month, if the timing uncertainties of the TOA of QPEs can be less than ∼ 100 − 150 seconds, the spinning of the MBH can be constrained. Orbitors with higher eccentricities, or QPE observations covering longer time span, can be used to set even stronger constraints on the spinning of MBH.
If the X-ray QPEs are driven by the star-disc collisions around the MBH, the time of arrival of each QPEs includes various general relativistic effects around the MBH. Our results suggest that long-time monitoring of X-ray QPEs can provide a unique probe of general relativistic effects, the spinning of MBH, tests of the no hair theorem, and also gravity theories.