The youngest of hot jupiters in action: episodic accretion outbursts in Gaia20eae.

Recent imaging observations with ALMA and other telescopes found widespread signatures of planet presence in protoplanetary discs at tens of au separations from their host stars. Here we point out that the presence of very massive planets at 0.1 au sized orbits can be deduced for protostars accreting gas at very high rates, when their discs display powerful Thermal Instability bursts. Earlier work showed that a massive planet modifies the nature of this instability, with outbursts triggered at the outer edge of the deep gap opened by the planet. We present simulations of this effect, finding two types of TI outbursts: downstream and upstream of the planet, which may or may not be causally connected. We apply our model to the outburst in Gaia20eae. We find that the agreement between the data and our disc thermal instability model is improved if there is a planet of 6 Jupiter masses orbiting the star at 0.062 au separation. Gaia20eae thus becomes the second episodically erupting star, after FU Ori, where the presence of a massive planet is strongly suspected. Future observations of similar systems will constrain the mode and the frequency of planet formation in such an early epoch.


INTRODUCTION
Ionisation of Hydrogen leads to a well known Thermal Instability (TI; Faulkner et al. 1983;Meyer & Meyer-Hofmeister 1984) in accretion discs, during which the inner disc switches between long periods of "quiescence" when the disc is neutral (cold), and shorter outbursts, when the disc is ionised (hot).With outbursts occurring multiple times a year in each of over 200 systems known, TI is the accepted theory of episodic accretion outbursts in dwarf novae (DNe; Dubus et al. 2018;Hameury 2020) and low mass X-ray binaries (LMXRBs; Coriat et al. 2012).These outbursts are unique in setting tight constraints on the usually very poorly known disc viscosity, indicating that α ∼ 0.01 (∼ 0.1) when H is the neutral (ionised) in the disc (Smak 1984;Lasota 2001).Recent 3D radiation and ideal magneto hydrodynamics (MHD) simulations of discs (e.g., Hirose 2015;Coleman et al. 2016;Scepi et al. 2018) reproduce such viscosity from first principles.
Cycles of low and high Ṁ are observed in accreting protostars as well, with the most powerful of them called FU Ori type (FUORs;Herbig 1977;Hartmann & Kenyon 1985).These outbursts are however disappointingly rare, with just a ⋆ sergei.nayakshin@le.ac.uk few known pre-2000.All of these outbursts also appeared unhelpfully long, e.g., ≳ 10 2 years (Hartmann & Kenyon 1996), preventing us from knowing how they end.The TI theory, re-scaled from stellar binaries to YSO, matches many of the observed "classic FUOR" characteristics (Clarke et al. 1990;Bell & Lin 1994;Clarke et al. 2005) but only if α is two orders of magnitude lower than in DNe/LMXRBs.While lower viscosity is plausible in quiescent cold discs due to non-ideal MHD effects (e.g., Gammie 1996;Gressel et al. 2015), this does not appear physically justifiable in outburst, prompting many to question the relevance of TI scenario for FUORs (e.g., Zhu et al. 2009;Armitage 2015).
Recent wide sky surveys uncovered a number of outbursts FUOR-like in luminosity but only ∼ 1 to a few years in duration.Nayakshin et al. (2024) (N24 hereafter) showed that these short outbursts are naturally produced by TI with the same viscosity as in DNe/LMXRBs, and presented statistical arguments that TI operates in the whole episodic accretors sample.However, N24 noted that additional physics must be present to (a) explain the long classic FUORs, and (b) account for suppression of most short TI outbursts since VVV/VVVX survey shows them far less frequent than theory predicts (Contreras Peña et al. 2024).This additional physics could be other disc instabilities (e.g., Armitage et al. 2001), stellar mass perturbers (e.g., Vorobyov et al. 2021), or massive gas planets (Vorobyov & Basu 2006, VB06 hereafter).
Studying TI in conjunction with other physics of YSO discs is natural (e.g., Pavlyuchenkov et al. 2023).TI operates in the inner ∼ 0.1 au.For example, in the Armitage et al. (2001) scenario an outburst is triggered at ∼ 1 au, prompting mass transfer towards the star at high rate (e.g., Bae et al. 2014;Martin & Lubow 2014).As the material propagates through the innermost 0.1 au, TI-outbursts are triggered (e.g., Appendix A1 in Nayakshin & Elbakyan 2024).If future work cements TI status as the first piece of the FUOR puzzle in place, then, following the lead of Scepi et al. (2019) in DNe, we can use TI to learn what else takes place in FUOR discs.
In this spirit, and motivated by the recently discovered short episodic accretor Gaia20eae (Cruz-Sáenz de Miera et al. 2022), here we study the model of Lodato & Clarke (2004) (LC04 hereafter) for FUORs.They showed that a gapopening (massive) planet migrating through the disc modifies TI outbursts significantly.While at separation a ≳ 1 au, the planet opens a practically impenetrable gap, and the disc downstream of the planet drains onto the star.The disc upstream of the planet, on the contrary, banks up and swells in mass.When the planet migrates to a ∼ 0.1 au, TI outbursts begin at the gap outer edge, and propagate outside-in rather than inside-out as in classic TI outbursts.LC04 used the same very low disc viscosity as Bell & Lin (1994) to obtain 100-year long outbursts1 .Here we use DNe-like viscosity values instead, and argue that effects described in LC04 are in action in Gaia20eae.

METHODS AND A FIDUCIAL CASE
Our 1D time-dependent planet-in-disc code was described in detail in Nayakshin et al. (2023) and in N24.Briefly, our numerical approach to evolving the disc surface density, Σ, is essentially identical to that of LC04, whereas our energy equation approach is more closely related to Bell & Lin (1994).We solve for Σ via where ν = αcsH is the kinematic viscosity (Shakura & Sunyaev 1973), cs and H are the disc midplane sound speed and vertical scale height, respectively; and Ω = (GM * /R 3 ) 1/2 is the Keplerian angular frequency.The last term in eq. 1 describes the angular momentum exchange between the disc and the planet via tidal torques; λ is specific torque on the disc material.The planet migration rate is set by the angular momentum conservation of the system.The equations for the evolution of the disc central temperature, Tc, are given in N24.Aside from disc viscosity, another distinction from LC04 is the disc heating due to dissipation of the tidal torque from the planet.This heating term, D tid (R), eq. 13 in Lodato et al. (2009) where a is planet-star separation, q = Mp/M * , ∆Rsm = max(H, RH ) is the distance within which both λ and D tid are smoothed to go to zero at ∆R = 0 (see LC04).One can show that D tid dominates over the viscous dissipation for where h = H/R.For q = Mp/M * = 5 × 10 −3 , eq. 3 yields ≈ 0.25 in outburst, but can be as large as 0.7 in quiescence.
As in N24 we compute disc spectra by integrating over the disc, assuming the local black-body emission spectrum with a time-dependent non power-law T eff (R), e.g., see Fig. 4c.We start with a planet of mass Mp = 5 MJ and initial separation a = 1 au inserted into a disc with cold and hot viscosity parameters αc = 0.02, α h = 0.1, and the mass feeding rate of the inner disc Ṁfeed = 7 × 10 −7 M⊙ yr −1 , added in a narrow ring close to the outer boundary of the grid.We take M * = 1.15M⊙ for Gaia20eae (Cruz-Sáenz de Miera et al.

2022).
Even though our α values are two orders of magnitude larger than that of LC04, our models reproduce most of their results.Note that since Ṁ = 3πνΣ in an unperturbed planetfree disc at R ≫ R * , our disc is ∼ 100 times lighter at the same Ṁ .Thus the planet completely dominates system evolution in the inner few au, even more so here than in LC04. to this and continues nearly classic planet-free TI outbursts, e.g., they are triggered at ∼ 2 stellar radii and propagate outward (see Fig. 4 in N24).In Fig. 1, snapshots at t = 1040 and 1059 show Σ just before and at the very end of an outburst.However, with time the inner disc drains onto the star.By t ≈ 2200 years, when a ≈ 0.4 au, the downstream Σ never reaches Σcrit, and TI outbursts stop.The planet is then slowly pushed towards the star while the downstream disc is completely emptied so an inner hole is formed (e.g., see t = 4, 000 in Fig. 1a).
The disc upstream of the planet, on the contrary, banks up and swells in mass.As in LC04, when the planet migrates to the inner ∼ 0.1 au, the peak in Σ just upstream of the planet edges very close to Σcrit, and the outbursts can restart there.The outbursts are triggered upstream of the planet and propagate in.The last snapshot (dashed green) in Fig. 1 shows Σ at the end of such an outburst.
Since we include D tid here, and due to a much higher α than in LC04, we find additional observable signatures of upstream TI outbursts.Fig. 2 shows mass fluxes in two important locations in the disc: Ṁ in the accretion rate on the star, while Ṁgap = −2πΣ(Ra)vR(Ra) is the mass flux through Ra, the radial bin containing the planet at that moment.Here vR is the radial velocity of the gas in that bin.These mass fluxes are shown at three moments illuminating the planet-disc system evolution best.The central legend in each panel shows the respective planet location.
In the top panel in Fig. 2, soon after the planet has been inserted into the disc, Ṁ is characteristic of the classic TI outbursts, whereas Ṁgap plunges to negligible levels as the gap is being excavated.As explained above, this cuts off the innermost disc of further fuel supply, and the downstream outbursts stop when the innermost disc runs out of matter; Ṁ then plunges to very small levels too.
This goes on until the planet is at a ≈ 0.08a au, when the disc just upstream of the planet is hot enough to experience local "upstream bursts", discussed more fully in §3.In these bursts the disc at the outer edge of the gap can heat up to T ≳ 10 4 K (Fig. 4b), and the gap region becomes leaky.The material then rushes through the gap into the downstream disc (cf. the middle panel in Fig. 2), but initially no TI activity takes place there because Σ < Σcrit.After a succession of "gulps" through which the downstream disc receives its mass from upstream of the planet, the downstream bursts return at t ≈ 4500 years.These bursts are different from the classic planet-free ones, though, because the mass supply into the inner disc is still controlled by the upstream bursts.We believe that LC04 did not see the gulps that fail to propagate a TI burst all the way to the star because their bursts were dozens of times longer, each bringing much more matter downstream from the planet than ours.
As the planet migrates further in, the upstream bursts become more and more powerful, and are eventually able to push through enough mass in each gulp for the downstream disc to experience a TI burst every time.The upstream and downstream outbursts are then synchronised in the sense that the former always causes the latter (bottom panel in Fig. 2).
While investigating parameter space of this scenario we came across outbursts with lightcurves resembling that of Gaia20eae.These outbursts start upstream of the planet, and the planet tends to be at a < 0.1 au.

THE MODEL FOR GAIA20EAE
To save computational time, our simulations aimed at addressing Gaia20eae start with the planet injected into the disc at a = 0.2 au at t = 0. To approximate the complete inner hole that the planet opens in more realistic simulations where it migrates from much further away, we set disc Σ to zero within 0.25 au.The planet is then pushed towards the star by the planet-disc interactions.By the time the planet is at a = 0.1 au the system adjusts to a structure resembling the one we found in longer time simulations discussed in §2.
We experimented with several free parameters of our model, i.e., Mp, viscosity parameters αc and α h , and Ṁfeed to find a combination of parameters that reproduce Gaia20eae outburst as closely as possible.Here we focus only on the "best case scenario" model, for which Mp = 6 MJ, αc = 0.04, α h = 0.1, and Ṁfeed = 7 × 10 −7 M⊙ yr −1 .Fig. 3 shows the observed lightcurve of Gaia20eae in two bands, W1 and r ′ , with symbols, and two model lightcurves in the same bands.Fig. 3a shows the model from N24 (their Fig. 8), reproduced here for comparison.This model has no embedded planet.Fig. 3b is our best case scenario from this paper.The outburst analysed here occurrs when the planet reaches separation a = 0.062 au.
The presence of the planet improves the agreement between the simulations and the observations in several ways: (i) The disc is brighter in mid-IR before the outburst.This is due to the banking up of material upstream of the planet (e.g., Fig. 3 in LC04).Since local energy liberation is proportional to Σ, the disc is brighter upstream of the planet than without it.
(ii) There is a small, ∼ 1 mag, "precursor burst" seen in r ′ that begins at t ≈ 2020.0,several months before the main outburst.The burst ends with a short term drop in r ′ at t ≈ 2020.3,just before the rise towards the main burst.These effects are related to the gap region being bridged and are explained together with point (iv) below.
(iii) The initial rise to the maximum light in r ′ is less steep than in the classic TI case.Nayakshin & Elbakyan (2024) obtained similar conclusions for the lightcurve rises in the visual B and V bands for the model scenarios (their §5.1 and 5.2).Note that this contrasts with LC04 conclusions, which may reflect the large difference in viscosity parameters between the two studies.
(iv) There is a ∼ 1 mag step-like drop in r ′ at t ≈ 2021.This drop is caused by the gap in the disc re-emerging during the cooling phase of the burst.This curtails the spatial extent of the region bright in the r ′ band.
For a deeper discussion of points (ii) and (iv), Fig. 4 shows several snapshots of the main disc variables around the time of the precursor burst.As in Fig. 1, when Σ > Σcrit, a classic TI outburst should commence.The precursor burst in Fig. 3b is however triggered at a much smaller Σ, at R ≈ 0.075 au, cf. the cyan curve in Fig. 4b.This happens due to the strong extra heating that the disc receives from the planet via tidal dissipation D tid (eq.2), which has a local maximum at the trigger point (eq.3).At larger R, D tid ∝ ∆R −3 , so falls very rapidly.D tid also falls strongly within the gap where Σ is small.We found that the precursor bursts disappear in our model if we set D tid = 0.
The precursor bursts are very distinct from the classic TI ones.The latter are triggered very close to the star, and hence lead to an increase in T eff from ∼ 3, 000 K to as much as ∼ 10, 000 K (e.g., Fig. 4 in Nayakshin & Elbakyan 2024).This fuels ∼ 5 mag brightness increase of the disc in the optical bands.The precursor bursts are triggered on the outer edge of the gap, factor of ∼ 5 further out in the disc.Therefore, in precursor bursts the local effective temperature cannot increase much.The disc luminosity in the visual does not increase significantly while the ionisation front propagates inward due to the strong advective cooling (cf.Figs.4b and c).
During the main outburst the gap region remains somewhat suppressed in Σ and cooler than the disc just outside; nevertheless, the material is able to pass through it relatively easily.Fig. 5 depicts the mass flow rate Ṁ onto the star and through the gap region, Ṁgap, vs time.We observe that beginning of the precursor burst coincides material upstrem of the planet starting to flow through the gap region into the inner disc.The main outburst commences when the ionisation front reaches the star and Ṁ suddenly skyrockets.Fig. 2 also shows that the mass flow rate through the gap is comparable to Ṁ over much of the 2020.This confirms that the pre-burst phenomenon is crucial to powering the whole of the outburst; it only appears as a small (1 mag) step-up in the r ′ because the region is not as bright in the optical as the hotter regions downstream of the planet during the main burst.
When the disc starts to cool, the upstream region falls back into quiescence first.This triggers re-emergence of the gap and termination of the gas flow through the gap at t ≈ 2020.8.The inner disc starves of fuel, which necessitates a drop in Ṁ onto the star seen in the black curve in Fig. 5, and a ∼ 1 mag drop in r' band shortly thereafter.The inner disc outburst continues until ≈ 2021.5 when it finally falls into quiescence.

DISCUSSION AND CONCLUSIONS
Here we confirmed results of Lodato & Clarke (2004) for TI outbursts in discs with massive gap-opening planets, and found new specific signatures of this scenario.Such planets isolate the inner disc from the mass reservoir -the outer disc -so the classic inside-out TI outbursts disappear together with the inner disc.When the planet migrates inside R ∼ 0.1 au, the outbursts restart, but are triggered upstream of the planet and propagate outside-in.
We used a much larger disc viscosity than did LC04, as is consistent with TI in stellar binary systems.Due to this, TI outbursts starting on the gap outer edge are much shorter and transfer much less mass into the inner disc.The mass of the latter hence increases in small "gulps" (cf.Fig. 2).A number of these is needed to refill the inner disc so that TI can restart there.Eventually the upstream and downstream discs come into equilibrium in which every gulp tips the inner disc over the instability threshold, and TI propagates all the way to the star.The outbursts in this case start with a small (in terms of the source magnitude change, say ∆r ′ ) step-up burst, which corresponds to the "gulp".This is then followed by the downstream burst with a much larger ∆r ′ .
We applied this scenario to the outburst of Gaia20eae (Cruz-Sáenz de Miera et al. 2022), finding a somewhat better fit to its lightcurve than without a planet.We note that the fit is far from perfect, but it is probably not realistic to expect anything better with a 1D modelling that employs a simple αprescription for viscosity.Further, some of the disagreements may be expected given missing physics in the model.For example, Fig. 9 in Bell (1999) shows how reprocessing of the radiation in the R ≳ 1 AU disc makes the disc significantly brighter during FUOR outbursts in IR; we do not include these effects.It is possible that disc irradiation could make our model brighter in W1 band, and make it more consistent with the observations of Gaia20eae.
Additionally, the stellar magnetosphere is neglected in the model (since magnetospheres appear to be crushed by the disc during FUOR outbursts, e.g., Hartmann & Kenyon 1996).However, closer to the end of the outburst the magnetosphere may recover.Preliminary numerical experiments show that this may explain the smooth decay into quiescence in Gaia20eae compared with theoretical model abrupt cutoff in Fig. 3.
Notwithstanding these model shortcomings, our calculations show that the existence of the gap opened by the planet in the disc has specific lightcurve signatures that may help to uncover massive planet presence in FUOR discs.Finally, another possible signature of planet presence is a short time scale quasi-periodic variability in the disc line spectra and  Step-up trigger The planet is at a=0.062 au   3. Before the outburst, a deep gap cuts the inner disc off from a steady mass supply from outside.The outburst is triggered just upstream of the planet (cyan curve).As the gas heats up, the material is able to cross the gap and fill the inner disc sufficiently for the outburst to reach all the way to the star.photometry, as is observed in FU Ori (e.g., Herbig et al. 2003;Clarke & Armitage 2003;Powell et al. 2012;Siwak et al. 2018).For a planet without mass loss, as considered here, the pumping mechanism of such variability would be planet-disc tidal interactions, especially if the planet gains eccentricity.Future 2D simulations are needed to evaluate such effects better.

AKNOWLEDGEMENT
We thank Cathie Clarke and Giuseppe Lodato for illuminating discussions of the problem.SN acknowledges the funding from the UK Science and Technologies Facilities Council, grant No. ST/S000453/1.This work was also supported by the NKFIH excellence grant TKP2021-NKTA-64.FCSM received financial support from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme (ERC Starting Grant "Chemtrip", grant agreement No 949278).

DATA AVAILABILITY
The data obtained in our simulations can be made available on reasonable request to the corresponding author.

Figure 1 .
Figure1.Disc Σ(R) for a selection of time moments in the fiducial case.Vertical lines mark the respective position of the planet.After the planet is introduced into the disc, it starves the disc downstream of matter, opening a complete inner hole.No TI outbursts occur in that configuration.However, when the planet migrates to ∼ 0.08 au, TI outburst restart, but this time they are triggered upstream of the planet.

Fig. 1 Figure 2 .
Figure2.Accretion rate onto the star (black) and through the gap region (red) vs time at three moments in the fiducial simulation.The large peaks in Ṁ and Ṁgap are TI bursts downstream and upstream of the planet, respectively.The planet gap region is impenetrable barrier in the top panel but becomes leaky when the planet migrates sufficiently in (the two bottom panels).Note that an upstream burst initiates each of the downstream bursts in the bottom panel.The time span in the middle panel is far longer than in the two other panels.See §2 for detail.

Figure 3 .
Figure 3. Observed photometry in two wavebands (symbols; WISE W1 and Sloan r' bands) and model (lines) for the outburst in Gaia20eae.Panel (a): the classic TI model from Fig. 8 in Nayakshin et al (2024); Panel (b): same but with an embedded planet, this paper.Two moments in time are shown, with the respective planet locations listed in the legend.The banking up of material upstream of the planet explains the gradual rise in W1, and the pre-burst in r' at early 2020.The drops in the lightcurve in the r' band seen at 2020.3 and 2021.0 correspond to the times when the flow bridges the gap, and when it re-opens, respectively.

Figure 4 .
Figure 4. Disc evolution trhough the early phases of the outburst shown in the right panel of Fig.3.Before the outburst, a deep gap cuts the inner disc off from a steady mass supply from outside.The outburst is triggered just upstream of the planet (cyan curve).As the gas heats up, the material is able to cross the gap and fill the inner disc sufficiently for the outburst to reach all the way to the star.

Figure 5 .
Figure5.Mass accretion rate onto the star (black) and through the planet location (red).The step-up 1 mag burst seen in Fig.3at t ≈ 2020 begins when the gap is bridged by the material crossing it from the exterior to the interior disc.As the disc cools, the gap closes, which induces a step down in Ṁ at around t ≈ 2021.