A Mildly Relativistic Outflow Launched Two Years after Disruption in the Tidal Disruption Event AT2018hyz

We present late-time radio/millimeter (as well as optical/UV and X-ray) detections of the tidal disruption event (TDE) AT2018hyz, spanning $970 - 1300$ d after optical discovery. In conjunction with earlier deeper limits, including at $\approx 700$ d, our observations reveal rapidly rising emission at $0.8-240$ GHz, steeper than $F_\nu\propto t^5$ relative to the time of optical discovery. Such a steep rise cannot be explained in any reasonable scenario of an outflow launched at the time of disruption (e.g., off-axis jet, sudden increase in the ambient density), and instead points to a delayed launch. Our multi-frequency data allow us to directly determine the radius and energy of the radio-emitting outflow, showing that it was launched $\approx 750$ d after optical discovery. The outflow velocity is mildly relativistic, with $\beta\approx 0.25$ and $\approx 0.6$ for a spherical and a $10^\circ$ jet geometry, respectively, and the minimum kinetic energy is $E_K\approx 5.8\times 10^{49}$ and $\approx 6.3\times 10^{49}$ erg, respectively. This is the first definitive evidence for the production of a delayed mildly-relativistic outflow in a TDE; a comparison to the recently-published radio light curve of ASASSN-15oi suggests that the final re-brightening observed in that event (at a single frequency and time) may be due to a similar outflow with a comparable velocity and energy. Finally, we note that the energy and velocity of the delayed outflow in AT2018hyz are intermediate between those of past non-relativistic TDEs (e.g., ASASSN-14li, AT2019dsg) and the relativistic TDE Sw\,J1644+57. We suggest that such delayed outflows may be common in TDEs.


INTRODUCTION
A tidal disruption event (TDE) occurs when a star wanders sufficiently close to a supermassive black hole (SMBH) to be torn apart by tidal forces, leading to the eventual formation of a transitory accretion flow (Rees 1988;Komossa 2015). Optical/UV and X-ray observations of TDEs are generally thought to track the mass fallback and accretion (e.g., Stone et al. 2013;Guillochon & Ramirez-Ruiz 2013). Radio observations, on the other hand, can reveal and characterize outflows from TDEs (Alexander et al. 2020), including the presence of relativistic jets (Zauderer et al. 2011;Giannios & Metzger 2011;De Colle et al. 2012).
Recently, two TDEs have been reported to show radio emission with a delay relative to the time of optical discovery. ASASSN-15oi was first detected ≈ 180 days after optical discovery with a luminosity that exceeded earlier radio limits (at 8, 23, and 90 days) by a factor of ≈ 20 (Horesh et al. 2021a). The radio emission subsequently declined until about 550 days, and then exhibited a second rapid rise with a detection at 1400 days with an even higher luminosity than the first peak; see Figure 2. iPTF16fnl was first detected ≈ 150 days after optical discovery, with a luminosity about a factor of 8 times larger than earlier limits (extending to 63 days) and appeared to slowly brighten to about 417 days (Horesh et al. 2021b). The initial abrupt rise in ASASSN-15oi seems distinct from the radio light curve of AT2019dsg, although both reach their peak radio luminosity on a similar timescale and at a similar level. The gradual rise and much lower peak luminosity of iPTF16fnl (≈ 10 37 erg s −1 ), on the other hand, may indicate that it is simply a less energetic example of typical radio-emitting TDEs.
Delayed radio emission from TDEs has been speculated to result from several possible effects. First, it may be due to a decelerating off-axis relativistic jet launched at the time of disruption (e.g., Giannios & Metzger 2011;Mimica et al. 2015;Generozov et al. 2017). Second, it may be due to an initial propagation of the outflow in a low density medium, followed by interaction with a significant density enhancement (e.g., Nakar & Granot 2007). Finally, it may be due to a delayed launch of the outflow compared to the time of disruption and optical/UV emission, for instance as might result from a state transition in the accretion disk (e.g., Tchekhovskoy et al. 2014) or the delayed accumulation of magnetic flux onto the black hole (e.g., Kelley et al. 2014). These scenarios can in principle be distinguished through a combination of detailed temporal and spectral information, that can be used to infer outflow properties such as radius, velocity, and energy, as well as the ambient density and its radial structure. In the case of ASASSN-15oi it has been speculated that a delayed outflow may best explain the radio data (Horesh et al. 2021a).
Against this backdrop, here we report the detection of rapidly rising radio emission from AT2018hyz (z = 0.0457) starting about 970 days after optical discovery, with a factor of 30 increase in luminosity compared to upper limits at about 700 days. Our extensive multi-frequency data, spanning 0.8 to 240 GHz, as well as optical/UV and X-rays, allow us to characterize the rising phase of emission in detail for the first time and hence to distinguish between possible scenarios for such a late rise. We find that the radio emission requires an energetic (E K ≈ 10 50 erg) mildly relativistic (β ≈ 0.2 − 0.6) outflow launched with a significant delay of about 700 − 750 days after optical discovery. An off-axis relativistic jet launched near the time of optical discovery, or a large density enhancement can be ruled out.
The paper is structured as follows. In §2 we describe our new late-time radio, mm, UV/optical, and X-ray observations, and in §3 we contrast the radio emission from AT2018hyz with those of previous TDEs. In §4 we model the radio spectral energy distribution and carry out an equipartition analysis to derive the physical properties of the outflow and environment. In §5 we describe the results for a spherical and a collimated outflow geometry. We discuss the implications of a delayed mildly-relativistic outflow in §6 and summarize our findings in §7.

Radio Observations
AT2018hyz was observed in targeted observations with AMI-LA 32 days after optical discovery leading to an upper limit of F ν (15.5 GHz) 85 µJy (Horesh et al. 2018), and with ALMA at 45 and 66 days to limits of F ν (100 GHz) 38 and 48µJy, respectively (Gomez et al. 2020). The location of AT2018hyz was also covered in the ASKAP Variables and Slow Transients Survey (VAST; Murphy et al. 2021) at 697 days with a limit of F ν (0.9 GHz) 0.9 mJy, and in the VLA Sky Survey epoch 2.1 (VLASS; Lacy et al. 2020) at 705 days with a limit of F ν (3 GHz) 0.45 mJy. We determine the flux density limits in the survey images with the imtool fitsrc command within the pwkit package 1 (Williams et al. 2017).
As part of a broader study of late-time radio emission from TDEs, we observed AT2018hyz with the Karl G. Jansky Very Large Array (VLA) at 972 days in C-band (Program ID 21A-303, PI: Hajela) and detected a source with a F ν (5 GHz) = 1.39 ± 0.02 mJy and F ν (7 GHz) = 1.00 ± 0.02 mJy (see Table 1). Following this initial detection we obtained multi-frequency observations spanning L-to K-band (≈ 1 − 23 GHz; Programs 21B-357, 21B-360 and 22A-458, PI: Cendes), which resulted in detections across the full frequency range. For all observations we used the primary calibrator 3C147, and the secondary calibrator J1024-0052. We processed the VLA data using standard data reduction procedures in the Common Astronomy Software Application package (CASA; McMullin et al. 2007), using tclean on the calibrated measurement set available in the NRAO archive. The observations and resulting flux density measurements are summarized in Table 1. Additionally, we use data collected by the commensal VLA Low-band Ionosphere and Transient Experiment (VLITE; Clarke et al. 2016) at 350 MHz during our multi-frequency observations (Table 1).
We also obtained observations with the MeerKAT radio telescope in UHF and L-band (0.8 − 2 GHz) on 2022 April 17 (1282 days; DDT-20220414-YC-01, PI: Cendes) and with the Australian Telescope Compact Array (ATCA) on 2022 May 1 (1296 days; Program C3472; PI: Cendes) at 2 − 20 GHz. For ATCA we reduced the data using the MIRIAD package. The calibrator 1934-638 was used to calibrate absolute flux-density and band-pass, while the calibrator 1038+064 was used to correct short term gain and phase changes. The invert, mfclean and restor tasks were used to make deconvolved wideband, natural weighted images in each frequency band. For MeerKAT, we used the flux calibrator 0408-6545 and the gain calibrator 3C237, and used the calibrated images obtained via the SARAO Science Data Processor (SDP) 2 . We confirmed via the secondary SDP products that the source fluxes in the MeerKAT images were ∼ 90% of the sources overlapping with the NRAO VLA Sky Survey (NVSS; Condon et al. 1998), and that frequency slices show a steady increase of flux within the MeerKAT frequency range (with a spectral index α ≈ −1).

Millimeter Observations
Following the initial VLA radio detection we also observed AT2018hyz with the Atacama Large Millimeter/submillimeter Array (ALMA) at 1141, 1201, and 1253 days (Project 2021.1.01210.T, PI: Alexander). These observations roughly coincide with the three multi-frequency VLA observations. The first and third ALMA observations were in band 3 (mean frequency of 97.5 GHz) and band 6 (240 GHz); the second epoch was in band 3 only. For the ALMA observations, we used the standard NRAO pipeline (version: 2021.2.0.128) in CASA (version: 6.2.1.7) to calibrate and image the data. We detect AT2018hyz in all observations with a rising flux density (Table 1). We note that in the first observation, the flux density is at least 30 times higher than the early upper limits.

X-ray Observations
We obtained a Director's Discretionary Time observation of AT2018hyz on 2022 March 19 (δt = 1253 days) with ACIS-S onboard the Chandra X-ray Observatory with an exposure time of 14.9 ks (Program 23708833, PI: Cendes). We processed the data with chandra_repro within CIAO 4.14 using the latest calibration files. An X-ray source is detected with wavdetect at the position of AT2018hyz with a net count-rate of (8.9 ± 2.5) × 10 −4 c s −1 (0.5 − 8 keV), with statistical confidence of 6 σ (Gaussian equivalent).
We extracted a spectrum of the source with specextract using a 1.5 radius source region and a source-free background region of 22 radius. We fit the spectrum with an absorbed simple power-law model (tbabs*ztbabs*pow within Xspec). The Galactic neutral hydrogen column density in the direction of AT2018hyz is N H,MW = 2.67 × 10 20 cm −2 (Kalberla et al. 2005), and we find no evidence for additional intrinsic absorption, to a 3σ upper limit of N H,int 2.8 × 10 22 cm −2 . The photon index is Γ = 1.5 ± 0.7 (1 σ confidence level), and the 0.3 − 10 keV unabsorbed flux is F x = 1.78 +0.35 −0.28 × 10 −14 erg s −1 cm −2 . The flux density at 1 keV is 4.03 ± 1.20 × 10 −6 mJy. The X-ray flux is a factor of about 2 lower than the flux measured with XRT in the first 86 days, F x = 4.1 +0.6 −0.4 × 10 −14 erg s −1 cm −2 (Gomez et al. 2020). This indicates that AT2018hyz has faded in X-rays over this time scale, and we discuss the implications of this in §4.3.

UV/Optical Observations
We obtained UV observations of AT 2018hyz from the UV/Optical Telescope (UVOT; Roming et al. 2005) on board the Neil Gehrels Swift observatory (Swift; Gehrels et al. 2004) on 2022 January 7 (1182 days post-discovery) to continue tracking the evolution of the UV light curve originally presented in Gomez et al. (2020). We measure an AB magnitude of m UVW2 = 20.20 ± 0.08 mag using a 5 aperture centered on the position of AT 2018hyz using the HEAsoft uvotsource function (Heasarc 2014). We correct this magnitude for Galactic extinction using the Barbary (2016) implementation of the Cardelli et al. (1989) extinction law (0.28 mag) and subtract the host contribution of 21.76 mag (Gomez et al. 2020) to obtain a final value of m UVW2 = 20.14 ± 0.08 mag.
We compare this measurement to the MOSFiT model of AT 2018hyz presented in Gomez et al. (2020), which fit all optical/UV data to about 300 days post discovery. The model predicts a UVW2 magnitude of 23.3 ± 0.2 at the time of our new observation, a factor of 18 times dimmer than observed. This points to excess emission in the UV band. . While the source is rising in all frequencies during the first radio detections, we find the source has begun to fade in L-band (1.4 GHz, yellow) and S-band (3.0 GHz, green) after ∼ 1250 days. In contrast, at higher frequencies such as C-band (5.5 GHz, light blue), X-band (9 GHz, dark blue), Ka-band (14 GHz, purple), pink), and in the millimeter band (97.5 GHz, brown, and 240 GHz, black) the source is still rising as roughly Fν ∝ t 5 through 1300 days. In the UHF-band (0.88 GHz, red) we see the source has risen in luminosity ∼ 2.25× from 1000 to 1280 days, but do not have enough sampling to establish whether it is decreasing.
Similarly, we measure the late time r-band magnitude of AT 2018hyz at ∼ 1030 days post-discovery by downloading raw ZTF images from the NASA/IPAC Infrared Science Archive 3 and combining all r-band images taken within ±12 days of MJD = 59501. We subtract a template archival pre-explosion ZTF image from the combined image using HOTPANTS (Becker 2015), and perform PSF photometry on the residual image. For calibration, we estimate the zeropoint by measuring the magnitudes of field stars and comparing to photometric AB magnitudes from the PS1/3π catalog. Corrected for galactic extinction, we find r = 21.83 ± 0.36 mag. At this phase, the MOSFiT model predicts a magnitude of r = 24.65 ± 0.30, a factor of 13 times dimmer than measured. We do not detect emission in g-band to a 3σ limit of 21.76 mag. This implies a late time color of g − r −0.1, redder than the latest measurements from Gomez et al. (2020) of g − r = −0.58 ± 0.23. 0.19 ± 0.17 1126 4.98 ± 0.05 9.18 ± 0.01 2.35 ± 0.03 0.06 ± 0.02 1199 8.22 ± 0.16 9.12 ± 0.03 2.27 ± 0.06 0.22 ± 0.07 1251 8.82 ± 0.10 9.17 ± 0.02 2.09 ± 0.03 0.14 ± 0.03 1282 8.83 ± 0.68 9.47 ± 0.04 2.24 ± 0.17 0.74 ± 0.14 Note-a The values for the first epoch are for a fixed value of p = 2.3, the average of the first two full SEDs at 1126 and 1199 days.
The radio light curves of AT2018hyz at frequencies of ≈ 0.9 − 240 GHz are shown in Figure 1. At all well-sampled frequencies we find a rapid rise at 950 days. At S-band (3 GHz) the flux density rises by at least a factor of 16 from the non-detection at 705 days to a peak at about 1250 days. This corresponds to a steep power law rise (F ν ∝ t α ) with α 4.8. Similarly, at C-band (5 − 7 GHz) we find a steady rise from about 1.4 mJy (972 days) to 7.8 mJy (1296 days) corresponding to α ≈ 6. A similarly steep rise is observed up to 240 GHz. Such a steep rise occurring across a large spectral range is not expected in any model of delayed emission due to an off-axis viewing angle, a decelerating outflow, or a rapid increase in the ambient density (e.g., Nakar & Piran 2011; see §5). Instead, the inferred steep power law rise indicates that the launch time of the outflow actually occurred much later than the time of optical discovery; for example, to achieve a power law rise of t 3 , as expected for a decelerating outflow in a uniform density medium, requires a delay launch of ∼ 600 days after optical discovery.
We note that at frequencies of 3 GHz, our latest observation indicates divergent behavior relative to the higher frequencies, with a pronounced decline in the flux density. For example, in L-band (1.4 GHz) we find a rapid decline from 8.7 to 5.3 mJy in the span of only 31 days (1251 to 1282 days). This differential behavior is due to rapid evolution in the shape of the spectral energy distribution (see §4.2).
In Figure 2 we show the radio light curve of AT2018hyz in the context of previous radio-emitting TDEs. The radio luminosity of AT2018hyz rapidly increases from 7 × 10 37 erg s −1 at ≈ 700 days to ≈ 2 × 10 39 erg s −1 at ≈ 1300 days, making it more luminous than any previous non-relativistic TDE. The rapid rise in AT2018hyz is even steeper than the second rising phase of ASASSN-15oi (see Figure 2; Horesh et al. 2021a), although the light curve of the latter contains only 2 data points (at 550 and 1400 days) and its actual rise may be steeper and comparable to AT2018hyz. We also note that due to the wide gap in radio coverage of AT2018hyz between about 80 and 700 days, as well as the relatively shallower early radio limits compared to ASASSN-15oi, it is possible to "hide" an initial bump in the light curve as seen in ASASSN-15oi at ≈ 180 − 550 days ( Figure 2); indeed, it is even possible that AT2018hyz had early radio emission comparable to that of AT2019dsg (Cendes et al. 2021a; Figure 2), which had a nearly identical radio peak luminosity and timescale to ASASSN-15oi, but a more gradual and earlier rise.
Finally, we note that the radio emission from AT2018hyz is still about a factor of 20 times dimmer than that of Sw J1644+57 at a comparable timescale (1300 days) and about 80 times dimmer than its peak luminosity ( Figure 2). Since the powerful outflow in Sw J1644+57, with an energy of ≈ 10 52 erg became non-relativistic at ≈ 700 days (Eftekhari et al. 2018), this again argues against an off-axis jet interpretation for the less luminous (and hence less energetic) radio emission in AT2018hyz; namely, in such a scenario the radio emission would have peaked significantly earlier and with a much higher luminosity.
In the subsequent sections we model the radio SEDs to extract the physical properties of the outflow and ambient medium, as well as their time evolution, and show that these confirm our basic arguments for a delayed outflow.

Modeling of the Radio Spectral Energy Distributions
The radio/mm SEDs, shown in Figure 3, exhibit a power law shape with a turnover and peak at ≈ 1.5 GHz through 1251 days. At 1282 days, however, the peak of the SED shifts upwards to ≈ 3 GHz. A rapid shift to a higher peak frequency is unprecedented in radio observations of TDEs. The power law shape above the peak is characteristic of synchrotron emission.
We fit the SEDs with the model of Granot & Sari (2002), developed for synchrotron emission from gamma-ray burst (GRB) afterglows, and previously applied to the radio emission from other TDEs (e.g., Zauderer et al. 2011;Cendes et al. 2021b), using specifically the regime 4 of ν m ν a : Figure 3. Radio to millimeter spectral energy distributions ranging from 972 days to 1282 days. The grey lines are representative fits from our MCMC modeling ( §4.1). It is apparent that the first four epochs exhibit a roughly constant peak frequency (≈ 1.5 GHz) and a steadily rising peak flux density, while the final epoch exhibits a rapid shift to a higher peak frequency of ≈ 3 GHz. VLA data at 972 days are combined with VAST data at 1013 days. VLA data at 1126 days are combined with ALMA data at 1141 days. VLA data at 1199 days are combined with ALMA data at 1201 days. VLA data at 1251 days are combined with ALMA data at 1253 days. MeerKAT data at 1282 days is combined with ATCA data at 1296 days.
We determine the best fit parameters of the model -F ν (ν m ), ν a , and p -using the Python Markov Chain Monte Carlo (MCMC) module emcee (Foreman-Mackey et al. 2013), assuming a Gaussian likelihood where the data have a Gaussian distribution for the parameters F ν (ν m ) and ν a . For p we use a uniform prior of p = 2 − 3.5. We also include in the model a parameter that accounts for additional systematic uncertainty beyond the statistical uncertainty on the individual data points, σ, which is this a fractional error added to each data point. The posterior distributions are sampled using 100 MCMC chains, which were run for 3,000 steps, discarding the first 2,000 steps to ensure the samples have sufficiently converged by examining the sampler distribution. The resulting SED fits are shown in Figure 3 and provide a good fit to the data.
Our first observation at 972 days only includes 5 and 7 GHz, but clearly points to an optically thin spectrum with a peak at 5 GHz. We combine this observation with a VAST detection of the source at 1013 days with a flux of 1.3±0.03 mJy at .89 GHz (Horesh et al. 2022), indicating the peak is between these two frequencies. We make the reasonable assumption that the lack of evolution in ν p between 1126-1251 days indicates no serious change at 972 days as well, and fix p = 2.30. We use these values to determine the physical properties of the outflow at 972 days.
From the SED fits we determine the peak frequency and flux density, ν p and F ν,p , respectively, which are used as input parameters for the determination of the outflow physical properties. The best-fit values and associated uncertainties are listed in Table 2. We find that ν p remains essentially constant at ≈ 1.3 − 1.6 GHz at 972 to 1251 days, while F ν,p increases steadily by a factor of 3.7. While the VLITE limits at 350 MHz lie above our SED model fits (Figure 8) we note that they require an SED peak at 0.5 GHz since otherwise these limits would be violated. In addition, the single power law shape above ν p indicates that the synchrotron cooling frequency is ν c 240 GHz at 1126 and 1251 days. For the SED at 1282 days we find that F ν,p has remained steady, while ν p increased by a factor of 2. We also note that the spectral index below the peak appears to be shallower than F ν ∝ ν 5/2 .
Given the unusual evolution to higher ν p in the latest epoch, we have also considered a model with two emission components, peaking at ν p,1 ≈ 1.5 and ν p,2 ≈ 3 GHz, in which the lower frequency component dominates the emission at early time, and the higher frequency component rises at later times. The results of this model are presented in the Appendix ( §7), but in the main paper we focus on the simpler single-component model.

Equipartition Analysis
Using the inferred values of ν p , F ν,p and p from §4.1, we can now derive the physical properties of the outflow and ambient medium using an equipartition analysis. In all epochs, we assume a mean value of p=2.3 in our calculations. We first assume the conservative case of a non-relativistic spherical outflow using the following expressions for the radius and kinetic energy ( where d L = 204 Mpc is the luminosity distance and z = 0.0457 is the redshift. The factors f A and f V are the area and volume filling factors, respectively, where in the case of a spherical outflow f A = 1 and f V = 4 3 × (1 − 0.9 3 ) ≈ 0.36 (i.e., we assume that the emitting region is a shell of thickness 0.1R eq ), while in the case of a jet 5 f A = (θ j Γ) 2 and f V = 0.27(θ j Γ) 2 . The factors of 4 1/13+2p and 4 11/13+2p in R eq and E eq , respectively, arise from corrections to the isotropic number of radiating electrons (N e,iso ) in the non-relativistic case. We further assume that the fraction of post-shock energy in relativistic electrons is e = 0.1, which leads to correction factors of ξ 1/13+2p and ξ 11/13+2p in R eq and E eq , respectively, with ξ = 1 + −1 e ≈ 11. We parameterize any deviation from equipartition with a correction factor = (11/6)( B / e ), where B is the fraction of post-shock energy in magnetic fields (here we use B = 0.1; see §4.3). Finally, χ e = (p − 2)/(p − 1) e (m p /m e ) ≈ 42.3, where m p and m e are the proton and electron masses, respectively.
Using R eq we can also determine additional parameters of the outflow and environment (Barniol Duran et al. 2013): the magnetic field strength (B), the number of radiating electrons (N e ), and the Lorentz factor of electrons radiating at ν a (γ a ): We note an additional factor of 4 and a correction factor of (γ a /γ m ) p−1 are added to N e for the non-relativistic regime (Barniol Duran, private communication); here we use γ m = max[χ e (Γ − 1), 2]. We determine the ambient density assuming a strong shock and an ideal monoatomic gas as n ext = N e /4V , where the factor of 4 is due to the shock jump conditions and V is the volume of the emitting region, V = f V πR 3 eq /Γ 4 , with f V as defined above.

Cooling Frequency and B
Using the inferred physical parameters we can also predict the location of the synchrotron cooling frequency, given by (Sari et al. 1998): Since ν c has a strong dependence on the magnetic field strength, measuring its location directly can determine B and whether the outflow deviates from equipartition. In Figure 4 we show our VLA+ALMA+Chandra SED at 1251 − 1253 days, along with our model SED from §4.1, which does not include a cooling break (dashed lines). This model clearly over-predicts the Chandra measurements. The steepening required by the Chandra data is indicative of a cooling break, which we model with an additional multiplicative term to Equation 1 of [1 + (ν/ν c ) s3(β3−β4) ] −1/s3 , where β 4 = −p/2 and we use s 3 = 10 ( Granot & Sari 2002;Cendes et al. 2021a). Fitting this model to the data we find ν c ≈ 10 13 Hz (solid lines in Figure 4). However, since the X-ray flux measured with Chandra is only a factor of 2 times fainter than the steady early-time X-ray flux, and hence can be due to a source other than the radio-emitting outflow, we consider it to be an upper limit on the contribution of the radio-emitting outflow. As a result, our estimate of ν c is actually an upper limit.
With the value of ν c determined, we adjust the value of B and solve Equation 7 after repeating the equipartition analysis (Equations 2 to 5) to account for the deviation from equipartition in those parameters. With this approach, , and models that include a cooling break to match the X-ray flux density (black; Equation 7). We find that νc ≈ 10 13 Hz.
we find that B ≈ 0.01. Given that the deviation is not significant, and that ν c is an upper limit (and hence B ≈ 0.01 is a lower limit) we conservatively assume equipartition ( B = 0.1) in our subsequent analysis and in Table 3. We emphasize that the change would be relatively minor if we adjusted B -at 1251 days, for the jetted model we find this would correspond with a radius decrease from log(R eq ) ≈ 18.21 to ≈ 18.16, and the energy would increase from log(E eq ) ≈ 49.80 to ≈ 49.98.

Spherical Outflow
We begin by investigating the properties of the outflow in the conservative case of spherical geometry. We summarize the inferred physical parameters for all epochs in Table 3. We find that the radius increases from log(R eq ) ≈ 17.22 to ≈ 17.57 between 970 and 1251 d, corresponding to a large velocity of β ≈ 0.28 over this time span. However, if we use the time of optical discovery as the outflow launch date, then the inferred velocity in the first epoch (972 d) is β ≈ 0.066 and in the fourth epoch (1251 d) it is β ≈ 0.11. This means that the assumption of a launch date that coincides with the optical discovery is incorrect. Instead, we find that the increase in radius during the first four epochs is roughly linear, and fitting such an evolution with the launch time (i.e., time at which R = 0) as a free parameter, we instead find t 0 ≈ 750 d; see Figure 5. Thus, the physical evolution of the radius during this period confirms our initial argument based on the rapid rise of the radio emission ( §3) that the outflow was launched with a substantial delay of about 2 years relative to the optical emission. The inferred outflow velocity using the four epochs is β ≈ 0.24, which is larger than the typical velocities of β ≈ 0.1 inferred for previous non-relativistic radio-emitting TDEs (Cendes et al. 2021a;Goodwin et al. 2022;Alexander et al. 2016Alexander et al. , 2017Anderson et al. 2019).
The outflow kinetic energy increases by a factor of about 5 from E K ≈ 1.1 × 10 49 to 5.8 × 10 49 erg between 970 and 1250 days. The rapid increase in energy indicates that we are observing the deceleration of the outflow. The kinetic energy is larger than that of previous radio-emitting TDEs by a factor of ≈ 4 (Goodwin et al. 2022;Anderson et al. 2019), although we note that this energy is a lower limit and would be higher if a deviation from equipartition is assumed.
Finally, we find that the inferred ambient density declines from n ext ≈ 1.9 to ≈ 1.0 cm −3 over the distance scale of 1.7 × 10 17 to 3.7 × 10 17 cm, or ρ(R) ∝ R −1 . This is similar to the rate seen in M87* and Sw J1644+57 and Sgr A*, and less steep than the density profiles inferred around previous thermal TDEs (see Figure 7 and citations therein).  (Table 3), assuming equipartition. We fit a linear trend (i.e., free expansion) to these data to determine the launch time of the outflow (green lines) and its uncertainty (grey shaded regions mark the 1σ range), excluding the observation at 1282 days (open circle) for reasons outlined in Section 4. We find that t 0,sphere = 750 +80 −127 days and t0,jet = 750 +73 −113 days relative to the time of optical discovery. We also mark the time of the final non-detection at 705 days for reference (vertical dashed line).
Combined with the mild decline in density with radius, this indicates that the late turn-on of the radio emission is inconsistent with a density jump.
To conclude, even in the conservative spherical scenario we find that radio emission requires a delayed, mildlyrelativistic outflow with a higher velocity and energy than in previous radio-emitting TDEs.

Collimated Outflow
In light of the large outflow velocity and kinetic energy in the spherical case, we also consider the results for a collimated outflow. In particular, we choose an outflow opening angle of 10 • , typical of GRB jets and the jet in Sw J1644+57 (Berger et al. 2012;Zauderer et al. 2013). For a collimated outflow the resulting radius (and hence velocity) are larger than in the spherical case, so we need to also solve for the Lorentz factor, Γ, which impacts the values of R eq (as well as the other physical parameters). We begin by using the launch time inferred in the spherical case, and then iteratively re-calculate Γ, R eq , and the launch date. In the process we also include the relevant modifications due to Γ in f A and f V , which also impact the value of R eq (Equation 2). We note that these corrections are relatively small but we account for them for completeness.
With this approach, we find that the launch date for a 10 • jet is t 0 ≈ 750 d post optical discovery. The resulting radii are a factor of 4.7 times larger than in the spherical case, leading to a mean velocity to 1251 days of β ≈ 0.57, or Γ ≈ 1.23. The kinetic energy is E K ≈ 6.3 × 10 49 erg at 1251 days. Finally, we find the density declines from n ext ≈ 1.0 to ≈ 0.5 cm −3 over the distance scale of 7.8 × 10 17 to 1.7 × 10 18 cm. This is less dense than the spherical case, but follows a similar profile of ρ(R) ∝ R −1 ( §5.1).

Off-Axis Jet and Other Scenarios
We can also consider the possibility that the late-time emission from AT2018hyz is caused by a relativistic jet with an initial off-axis viewing orientation. The radio emission from an off-axis relativistic jet will be suppressed at early times by relativistic beaming, but will eventually rise rapidly (as steep as t 3 ) when the jet decelerates and spreads.
The time and radius at which the radio emission will peak are given by (e.g., Nakar & Piran 2011): where β 0 is the initial velocity, which for an off-axis relativistic jet is β 0 = 1. Using the equipartition parameters in the spherical case, we find t dec ≈ 25 d and R dec ≈ 8.3 × 10 16 cm. The value of t dec is substantially smaller than the observed delay of ∼ 10 3 d (as is R dec ). This agrees with our general argument in §3 that the lower radio luminosity of AT2018hyz and its much later appearance compared to the radio emission of Sw J1644+57 (which became nonrelativistic at ≈ 700 d) argues against an off-axis jet launched at the time of disruption. Moreover, the radio emission from an off-axis jet is expected to rise no steeper than t 3 (Nakar & Piran 2011), whereas in AT2018hyz the emission rises as t 5 if the outflow was launched at the time of disruption. We can consider the possibility of an off-axis jet expanding initially into a low-density cavity, followed by a denser region (thus delaying t dec to a longer timescale as observed). However, we can rule out this model because the observed rise in radio emission spans ≈ 300 d, while Equation 8 indicates deceleration over a timescale of t dec ≈ 30 d even if the time at which deceleration starts is itself delayed. Similarly, if there was initially a higher density environment which we did not capture in our observations, this would only correspond with a faster t dec .
We also consider the hypothesis that the rise in radio emission is caused by unbound material from the initial disruption colliding with a surrounding dense circumstellar material (CSM). Theoretical modeling of such unbound debris indicates the fastest speeds reached are v ≈ 0.03c (Yalinewich et al. 2019;Guillochon et al. 2016), which is significantly smaller than what we infer in both the spherical and jetted cases for AT2018hyz. We thus conclude this scenario is unlikely.
Finally, we considered a model in which the change in the SED properties during the latest epoch could be due to a combination of two outflows, with one dominating at a lower frequency and fading and the second dominating at higher frequencies and rising (see Appendix). Such a model may be expected if internal shocks within the outflow are leading to dissipation at more than one radius. However, we find that such a model requires the initial emission component to decline very rapidly (effectively turn off), while the later emission component has to rise by more than an order of magnitude in only 30 days between about 1250 and 1280 days. Such rapid evolution does not seem feasible, even in the context of AT2018hyz.
6. DISCUSSION AND COMPARISON TO PREVIOUS RADIO-EMITTING TDES

Outflow Kinetic Energy and Velocity
In Figure 6 we plot the kinetic energy and velocity (Γβ) of the delayed outflow in AT2018hyz in comparison to previous TDEs for which a similar analysis has been carried out, using the highest energy inferred in those sources (Cendes et al. 2021a;Stein et al. 2021;Alexander et al. 2016Alexander et al. , 2017Anderson et al. 2019;Goodwin et al. 2022;Zauderer et al. 2011;Cendes et al. 2021b). We find that at its peak, in the spherical case the energy is ≈ 4 times larger and the velocity is ≈ 2 times faster than in previous non-relativistic TDEs. If we compare to ASASSN-15oi specifically, using the observation with the highest peak frequency and peak flux (i.e., 182 days post optical discovery), with B = 0.1 and p = 2.4 (which best fits the SED; see Horesh et al. 2021a), we obtain an energy for ASASSN-15oi that is ≈ 23 times lower than in AT2018hyz. To infer the velocity in ASASSN-15oi, we subtract 90 days from the date of the observation (the last date of non-detection; see Horesh et al. 2021a) to find β ≈ 0.13, which is lower than β ≈ 0.25 for AT2018hyz.
If the outflow in AT2018hyz is collimated, then its velocity is significantly higher than in previous non-relativistic TDEs, placing it in an intermediate regime with the powerful jetted TDEs such as Sw J1644+57.

Circumnuclear Density
In Figure 7 we plot the ambient density as a function of radius (scaled by the Schwarzschild radius) for AT2018hyz and previous radio-emitting TDEs. Here we use M BH ≈ 5.2 × 10 6 M for AT2018hyz, as inferred by Gomez et al. (2020). We find that the density decreases with radius, and is consistent with the densities and circumnuclear density profiles of previous TDEs, including Sw J1644+57 Crucially, we do not infer an unusually high density, which might be expected if the radio emission was delayed due to rapid shift from low to high density.  Figure 6. The energy/velocity for AT2018 spherical case and jetted case; note in both cases the energy increases over time. We include non-relativistic TDEs assuming a spherical outflow (Cendes et al. 2021a;Stein et al. 2021;Alexander et al. 2016Alexander et al. , 2017Anderson et al. 2019;Goodwin et al. 2022), and Sw J1644+57 which launched a relativistic jet (Zauderer et al. 2011;Cendes et al. 2021b). In the case of AT2019dsg, we display the highest energy inferred in the system. In the case of AT2019azh, which showed significant fluctuations in the final energy due to changes in p in some later observations (Goodwin et al. 2022), we show the energy and velocity at the peak flux in the luminosity curve adjusted to B = 0.1. For ASASSN-15oi, we use the observation with the highest peak frequency and peak flux (182 days) B = 0.1 and p = 2.39, which best fit the SED at the highest peak flux observation at 182 days post-disruption, and infer the velocity by subtracting 90 days (the last date of non-detection; see Horesh et al. 2021a). We find AT2018hyz is more energetic than the non-relativistic outflow TDEs, and has a higher velocity.

Comparison to Other TDEs with Late Radio Emission
Two previous TDEs with delayed radio emission have been published recently. The radio emission from iPTF16fnl is detected on a much earlier timescale than in AT2018hyz and rises more gradually to a peak luminosity of only ≈ 10 37 erg s −1 (Horesh et al. 2021b). This is almost a factor of 100 less luminous than AT2018hyz. We do not consider the radio emission from iPTF16fnl to be similar in nature to AT2018hyz.
On the other hand, ASASSN-15oi exhibits two episodes of rapid brightening, at ≈ 200 and ≈ 1400 d (Horesh et al. 2021a). While the second brightening is not well characterized temporally or spectrally, it has a comparable rise rate and luminosity to AT2018hyz. We speculate that it may be due to a delayed outflow with similar properties to that of AT2018hyz, including a delay of several hundred days, which would make it distinct from the first peak in the ASASSN-15oi light curve. . The circumnuclear density profile derived from various TDEs including AT2018hyz, normalized to the Schwarzschild radius of the SMBH at each host galaxy's center. AT2018hyz's host galaxy is lower density, similar to that seen in the jetted TDE SwJ1644+57 (Berger et al. 2012;Zauderer et al. 2013;Eftekhari et al. 2018) and M87 (Russell et al. 2015). We also include non-relativistic TDEs (e.g. ASASSN-14li, Alexander et al. 2016, AT2019dsg Cendes et al. 2021a, AT2019azh (Goodwin et al. 2022, ASASSN-15oi (Horesh et al. 2021a)), and for the Milky Way (Baganoff et al. 2003;Gillessen et al. 2019). Other than for AT2019dsg and Sw J1644+57, where the cooling break is detected, we assume equipartition in all radio TDEs, so their densities are lower limits.
There are at least two broad possibilities for the origin of the delayed mildly relativistic outflow, both of which connect to its assumed origin in a fast disk wind or jet (hereafter, collectively referred to as jet) from the innermost regions of the black hole accretion flow.
One possibility is that the jet was weak or inactive at early times after the disruption and then suddenly became activated at ≈ 750 d. Such sudden activation could result from a state change in the SMBH accretion disk, such as a thin disk that transitioned to a hot accretion flow. This is predicted to occur -in analogy with models developed to explain state changes in X-ray binaries -when the mass accretion rate falls below a few percent of the Eddington rate (e.g., Tchekhovskoy et al. 2014). The optical light curve of AT2018hyz extends to ≈ 800 d (Gomez et al. 2020;Hammerstein et al. 2022;Short et al. 2020), overlapping the time at which we estimate the outflow was launched. When examining the data from ZTF at later times, we find there is an excess ( §2). Based on analytical modeling of the optical data using the MOSFiT code (Gomez et al. 2020) we estimate that the mass accretion rate at the time the radio outflow was launched is ≈ 0.05Ṁ Edd , consistent with the possibility of a state change.
An alternative explanation is that powering a relativistic jet via the Blandford-Znajek process requires a strong magnetic flux threading the black hole horizon. The original magnetic field of the disrupted star in a TDE is not expected to contain a strong enough magnetic field to power a relativistic jet (Giannios & Metzger 2011), which requires an alternative origin. The first possibility is that the magnetic flux could be generated through a dynamo by the accretion disk itself; Liska et al. (2020) found that it may take only ∼ 10 days to generate poloidal flux from the toroidal field through a dynamo effect once the disk is sufficiently thick, thereby connecting jet production to the disk getting thinner as the accretion rate drops. Alternatively, Tchekhovskoy et al. (2014) and Kelley et al. (2014) suggest that the required magnetic flux may originate from a pre-existing AGN disk, which is 'lassoed' in by the infalling fall-back debris; because the matter falling back at later and later times in a TDE reaches larger and larger apocenter radii, depending on the radial profile of the magnetic flux in the pre-existing AGN disk, this could delay the jet production.
It is also possible instead that the jet has been present for the entire duration of the TDE. However, due to the combination of the high density of the large cloud of circularizing TDE debris (e.g., Bonnerot et al. 2022) and the potential for jet precession (e.g., due to misalignment of the disk angular momentum relative to the black hole spin axis; e.g., Stone & Loeb 2012), the jet is initially choked. At later times, as the accretion rate and gas density surrounding the black hole drop, eventually the jet is able to propagate through the debris cloud and escape.
Planned additional monitoring of AT2018hyz may more clearly elucidate the mechanism responsible for the delayed launch.

CONCLUSIONS
We presented the discovery of late-and rapidly-rising radio/mm emission from AT2018hyz starting at about 970 d post optical discovery, and extending to at least 1280 days. The radio emission is more luminous than previous non-relativistic TDEs, but still an order of magnitude weaker than the relativistic TDE Sw J1644+57. The rapid rise in luminosity coupled with the slow spectral evolution to 1250 days point to a decelerating outflow. Basic modeling assuming energy equipartition indicates that the outflow was launched ≈ 700 − 750 d after optical discovery with a velocity of β ≈ 0.25 (spherical outflow) or up to ≈ 0.6 (10 • jet). The outflow kinetic energy is at least 6 × 10 49 erg. This is the first case of a delayed mildly-relativistic outflow in a TDE, and its energy is in excess of all previous non-relativistic TDEs. On the other hand, we find that the density of the circumnuclear environment is typical of previous TDEs, indicating that interaction with a dense medium is not the cause for the long delay. We similarly show that an off-axis jet cannot explain the late rising radio emission.
With planned continued observations of AT2018hyz we will monitor the on-going evolution of the outflow and of the circumnuclear medium. We note that the discovery of such late-time emission indicates that delayed outflows may be more common than previously expected in the TDE population. A systematic study of a much larger sample of TDEs will be presented in Cendes et al. in preparation. We thank Emil Polisensky for his assistance with VLITE data, and the VLA, MeerKAT, and Chandra observatory staff for prompt responses and observations for our DDT requests. In particular, we thank Drew Medlin for his assistance with the calibration and imaging of our wide-field VLA data. The Berger Time-Domain Group at Harvard is supported by NSF and NASA grants. This paper makes use of the following ALMA data: ADS/JAO.ALMA#2021.1.01210.T. TL acknowledges support from the Radboud Excellence Initiative. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada), MOST and ASIAA (Taiwan), and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ. The National Radio Astronomy Observatory is a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc. The scientific results reported in this article are based in part on observations made by the Chandra X-ray Observatory, which is operated by the Smithsonian Astrophysical Observatory for and on behalf of the National Aeronautics Space Administration under contract NAS8-03060. The MeerKAT telescope is operated by the South African Radio Astronomy Observatory, which is a facility of the National Research Foundation, an agency of the Department of Science and Innovation. The Australia Telescope Compact Array is part of the Australia Telescope National Facility (https://ror.org/05qajvd42) which is funded by the Australian Government for operation as a National Facility managed by CSIRO.  Figure 8. Radio spectral energy distributions for our VLA data as in Figure 3, but for a two outflow model. Here, Outflow 1 has a peak near ∼ 1 GHz (green lines), Outflow 2 has a peak near ∼ 3 GHz (orange lines), and they combine to give a total outflow (black lines).
As noted in §4, given the unusual nature of the spectral evolution between the latest two epochs, with the peak frequency increasing by a factor of 2 in about 31 days, we consider here a model with two outflows each generating its own synchrotron emission, and peaking at ≈ 1.5 GHz (Outflow 1) and ≈ 3 GHz (Outflow 2). In this scenario, Outflow 1 dominates the emission at ≈ 970 − 1250 days, while Outflow 2 becomes more dominant at ≈ 1280 days.
To fit the data we use the same procedure described in §4. We set the MCMC priors to ensure that the time evolution of ν p,1 and ν p,2 is to either remain constant or decline, and that the peak flux density of Outflow 2 increases with time.
We first model the first and final well sampled epochs (1126 and 1282 days) to constrain the range of ν a,1 and ν a,2 . We then model the intermediate timescales; we do not model the earliest epoch at 972 days since the SED is not well sampled enough for this two-component analysis. The results of the models are shown in Figure 8, and the parameters of the two outflows are provided in Table 4.
We find that Outflow 1 dominates the emission in the first three epochs, but has to rapidly fade from about 8 mJy to 3.6 mJy in the final epoch. Outflow 2 provides minimal contribution in the first three epochs, and then rises rapidly from about 0.4 mJy to 6.1 mJy in about 30 days. This would effectively correspond to Outflow 2 being generated at a later time, at about 1250 − 1280 days.