High-resolution VLBI Observations of and Modeling the Radio Emission from the Tidal Disruption Event AT2019dsg

When a star passes within the tidal radius of a galactic supermassive black hole (SMBH of mass M_• = 10⁶–10⁸ M_⊙, where M_⊙ is the solar mass), its destruction by tidal forces leads to the fallback accretion of stellar debris (Rees 1988; Evans & Kochanek 1989). This usually powers a thermal flare in the optical, ultraviolet (UV), and X-ray bands, with a luminosity comparable to bright supernovae at ≈10⁴³–10⁴⁴ erg s⁻¹ (e.g., Saxton et al. 2020). The optical/UV luminosity generally declines with a power-law index of −5/3 (e.g., Saxton et al. 2020) over months to years, indicative of fallback accretion (Evans & Kochanek 1989; Phinney 1989). However, there may be small or significant deviations from this trend (e.g., van Velzen et al. 2021). Further, the later-phase (months to years) decline in the X-ray and UV bands may be more shallow (with an index of −5/12) and is indicative of the prevalence of disk emission (Lodato & Rossi 2011). Outflows can be launched as a result of stream-stream collisions involving the interaction between the orbiting stellar debris streams and material accreting onto the SMBH (e.g., Jiang et al. 2016; Lu & Bonnerot 2020). These can also originate as radiation-pressure-driven disk winds (Strubbe & Quataert 2009). Upon interaction with the surrounding circumnuclear medium (CNM), the outflow produces a nonthermal, predominantly synchrotron afterglow emission (e.g., Roth et al. 2020; Alexander et al. 2016). Multiwavelength observations of the afterglow can help identify the onset of this interaction and constrain the energy output from the tidal disruption event (TDE) and the properties of the CNM (e.g., Alexander et al. 2020).

The transient AT2019dsg was discovered on 2019 April 9 by the Zwicky Transient Factory (Nordin et al. 2019) and identified as a TDE in a host galaxy at a redshift of z = 0.051 (Nicholl et al. 2019; van Velzen et al. 2021), making it a rarely seen nearby TDE. Stein et al. (2021) reported a putative association of AT2019dsg with the 0.2 peta electron volt (PeV) energy neutrino event IC191001A, which was detected by the IceCube Neutrino Observatory on 2019 October 1 (175 days after the TDE discovery; IceCube Collaboration 2019) and is expected to be of astrophysical origin. A multiwavelength analysis in the same study (Stein et al. 2021) indicates a highly luminous TDE (optical peak luminosity of ≈3.5 × 10⁴⁴ erg s⁻¹) and suggests the presence of an active central engine embedded in an optically thick photosphere, which powers subrelativistic outflows and provides a suitable site for the neutrino production. Optical polarimetric observations (Lee et al. 2020) show a decrease in the linear polarization degree from 9.6% near the optical peak around 2019 May 17 to 2.7% on 2019 June 20 during the declining phase; this variation is attributed to a nonisotropic accretion disk or a relativistic jet. A concordance scenario (Winter & Lunardini 2021) proposes the presence of a relativistic jet that participates in the acceleration of protons and provides the necessary site for neutrino production. Alternatively, an off-axis jet may host sites for hadronic interactions and neutrino production (Liu et al. 2020). An analysis of multiwavelength radio observations of AT2019dsg (Cendes et al. 2021a) suggests that the emission originates from the interaction of a nonrelativistic outflow with the surrounding medium. Based on the relatively low outflow velocity and energy output in the context of TDEs and Type Ib/c supernovae, the study of Cendes et al. (2021a) suggests that the neutrino association is less likely. Matsumoto et al. (2022) probe the nature of the radio emission through the associated synchrotron self-absorption, and use calculations based on an equipartition model to infer the presence of a freely expanding or decelerating nonrelativistic outflow. The complex observational signatures and associated inferences are thus inconclusive as to the origin of the nonthermal emission, which may originate from either a relativistic jet or a nonrelativistic outflow. Clearly distinguishing between them is crucial to constrain the TDE properties (physical and geometric properties of the black hole–disk system and the disrupted stellar system), the nuclear environment of the host galaxy, and the central engine activity and to clarify the neutrino production mechanisms.

Very long baseline interferometry (VLBI) observations at milliarcsecond (mas) resolutions can help measure the size of the radio afterglow, constrain the proper motion of the emitting component (e.g., Yang et al. 2016; Mattila et al. 2018; Mohan et al. 2020), and provide inputs to discriminate the TDE and CNM properties. Of those TDEs that are detected in radio (e.g., Alexander et al. 2020; Horesh et al. 2021), only a few have been imaged using VLBI, with large differences in radio luminosity evolution, environmental properties, and the level of collimation of the outflow. Besides, extensive monitoring campaigns have only been conducted on Swift J1644+5734 (e.g., Levan et al. 2011; Zauderer et al. 2011; Berger et al. 2012; Zauderer et al. 2013; Yang et al. 2016; Eftekhari et al. 2018; Cendes et al. 2021b) and Arp 299-B AT1 (e.g., Mattila et al. 2018) owing to the presence of a relativistic jet, ASASSN-15oi (Horesh et al. 2021) due to a peculiarly delayed radio flaring (over timescales of months to years), ASASSN-14li (e.g., Romero-Cañizales et al. 2016; van Velzen et al. 2016; Bright et al. 2018; Alexander et al. 2016) being relatively close in proximity (z = 0.02), and AT2019dsg owing to its persistent radio emission and potential neutrino association (Stein et al. 2021; Cannizzaro et al. 2021; Cendes et al. 2021a; Matsumoto et al. 2022). A multiwavelength observation campaign was carried out to study AT2019dsg and its host galaxy (Cannizzaro et al. 2021). The results include the inference of a central SMBH of 5.4 × 10⁶ M_⊙, derived from the optical spectroscopy. In addition, they report the detection of a compact source with no significant displacement (1σ astrometric uncertainty of ≈6 mas) from the phase center, as inferred from the 5 and 1.4 GHz radio monitoring observations made by the electronic Multi-Element Remotely Linked Interferometer Network (e-MERLIN) up to 180 days after the TDE.

We employ the VLBI technique to monitor any emission structure change of AT2019dsg at high resolutions, motivated by potentially identifying a (collimated) relativistic jet and determining the source size. This makes it one of only three thermal-emission-dominated TDEs with VLBI monitoring, the others being ASASSN-14li (Alexander et al. 2016) and CNSS J0019+00 (Anderson et al. 2020). Our observations were carried out with the European VLBI Network (EVN) at 5 GHz spanning three epochs at 200, 216, and 304 days after the TDE, covering the declining phase of the source radiative evolution. As the source was successfully detected in all three epochs, this provides unique constraint information: flux densities during the declining phase, and kinematics including source size and astrometry to infer any significant proper motion. The first and third observations were made in the standard disk-recording mode, involving a maximum of 20 telescopes. The second observation was made in the electronic-VLBI (e-EVN) mode. Details of the observation logs are presented in Table 1, and the data reduction is presented in Appendix A.

Figure 1 shows the EVN images of AT2019dsg during the three epochs and the constraints on the peak emission positions with respect to the phase center. The astrometric parameters involving the calibrators and target source including errors are presented in Table 2. With a maximum expansion distance of 0.32 ± 0.22 mas (1σ) that spans 123.7 days, the proper motion is 0.94 ± 0.65 mas yr⁻¹ (3.2 ± 2.2 c). ⁴ At the 1σ level of uncertainty, we cannot rule out a relativistic expansion. However, as the proper motion ≲1.4σ, it may be considered to be statistically insignificant (with a requirement to be ≥3σ). A discernible proper motion may be expected in the case of an off-axis jet. For an on-axis view, Doppler beaming of emission along the observer's line of sight is expected to result in high brightness temperatures (≥10¹⁰ K, the case for blazars and jetted active galactic nuclei, e.g., Cheng et al. 2020). Neither of these scenarios are applicable here owing to the absence of a statistically significant proper motion (large uncertainties at the 3σ level) and brightness temperatures of ∼10⁹ K and ≥10⁸ K during the first and third observation epochs (see Table 3). These and the detection of an unresolved compact structure indicate the presence of a nonthermal component, and are consistent with afterglow emission from the interaction of an outflow with the surrounding CNM, which produces an expanding forward shock.

Zoom In Zoom Out Reset image size

The compiled 5 GHz light curve consists of flux densities reported in Stein et al. (2021) and Cannizzaro et al. (2021), as well as the measurements from our EVN observations. A smooth broken power law is used to fit the light curve (see Appendix B and Table 4 for the fit-parameter estimates), presented in Figure 2. The peak flux density ${F}_{{\nu }_{p}}$ is inferred to be 1.18 ± 0.18 mJy at a time t_p of 152.8 ± 16.2 days. The power-law indices during the rising and declining phases are α₁ = 2.4 ± 0.9 and α₂ = −2.6 ± 1.3, respectively. The peak flux density, associated time, and the power-law temporal indices are used as inputs in the modeling and interpretation. For a spherically expanding forward shock with a constant electron-injection energy index p, a time dependence of size R ∝ t^α , a magnetic field strength B ∝ R⁻² ∝ t^−2α and the proportionality constant (for the number of electrons per unit of energy within the region) N₀ ∝ B² ∝ t^−4α result in optically thick and thin flux densities that scale as ∝ t^3α and ∝ t^−(2+p)α respectively. The expansion velocity v ∝ t^α−1; comparing the measured maximum speed of 0.57 c with that inferred from the VLBI measured size of 0.17 mas (0.17 pc) and the 5 GHz flux-density peak time t_p of 152.4 days, α ≥ 0.59, with the lower limit implying a deceleration of the expanding shock. This is consistent with the decelerating solution (index of 0.63) based on an equipartition analysis as inferred in Matsumoto et al. (2022) for a radio-emitting outflow interacting with the surrounding CNM. Using α = 0.59, and p = 2.7 ± 0.2 inferred in Cendes et al. (2021a) under the assumption of a spherically expanding outflow, the optically thin flux density scales as ∝ t^−2.8, consistent with the measured α₂ = −2.6 ± 0.3. Hence, we adopt these as physical parameters for subsequent calculations in the radiative models. The optically thin spectral index α_s depends on the observation frequency as ∝ ν^−(p−1)/2 ∝ ν^−0.85, indicating that the emission is dominated by a steep-spectrum component. The radio luminosity at 5 GHz is then $\nu {L}_{\nu }\approx 4\pi {D}_{L}^{2}{F}_{{\nu }_{p}}{(1+z)}^{-(1-{\alpha }_{s})}\approx 3.6\times {10}^{38}$ erg s⁻¹. This makes AT2019dsg relatively radio-quiet (e.g., Alexander et al. 2020) and comparable to ASASSN-14li (Alexander et al. 2016, 2020) and CNSS J0019+00 (Anderson et al. 2020), whose radio emission is thought to be associated with nonrelativistic outflows.

Zoom In Zoom Out Reset image size

Table 4. Fit Parameters for the Model ${F}_{\nu }={F}_{{\nu }_{0}}+{F}_{{\nu }_{1}}{({(t/{t}_{p})}^{5{\alpha }_{1}}+{(t/{t}_{p})}^{5{\alpha }_{2}})}^{-1/5}$

Parameter	Range	Fit Value
Normalization (constant) ${F}_{{\nu }_{0}}$ (mJy)	0.00–0.50	0.17 ± 0.17
Normalization (shape) ${F}_{{\nu }_{1}}$ (mJy)	0.00–3.00	1.18 ± 0.18
Rising slope α1	0.00–5.00	2.38 ± 0.93
Declining slope α2	−10.00–0.00	−2.56 ± 1.31
Peak time (days)	100.00–300.00	152.80 ± 16.19
Peak flux density ${F}_{{\nu }_{p}}$ (mJy)	⋯	1.19 ± 0.18

Note. The large uncertainties in the peak time and declining slope are a consequence of sparse data points in these evolution phases.

Download table as:  ASCIITypeset image

For afterglow emission from a relativistically expanding region, the expected time-dependent angular size is ${\theta }_{A}\propto {({E}_{j,52}/{n}_{0})}^{1/4}{t}^{1/4}$ and $\propto {({E}_{j,52}/{n}_{0})}^{1/2}{t}^{1/2}$ for a constant density and a stratified surrounding medium respectively (see Equations (C4) and (C6)), where E_j,52 is the total kinetic energy of the outflow scaled in units of 10⁵² erg and n₀ is the number density at a distance of 1 pc from the central engine. For the VLBI sessions 1 and 3, the constrained size θ_A,obs is ≤0.17 mas and ≤0.44 mas, respectively. In comparison with the above models, the best estimate gives (E_j,52/n₀) ≤ 0.08 for session 1 (200 days post TDE) and a consequent Lorentz factor Γ_sh ≤0.80 (see Equations (C3) and (C5)), indicative of a nonrelativistic expansion in AT2019dsg, consistent with the above expectation. The VLBI measured flux densities and constraint on source size thus play a complementary role in enabling inferences on the emission source nature, and help reduce the number of free parameters in the synchrotron radiation model. The lack of VLBI/high-resolution radio observations during the early epochs, and the consequent unavailability of astrometric accuracy necessitates the use of a model-based approach to infer the possibility of a relativistic expansion. This may last for ≤18 days (relativistic phase marked by a Lorentz factor Γ ≤ 2; see Equations (C3) and (C5)) based on the above expansion models; a nonrelativistic expansion is thus more likely at ≥18 days, and an appropriate radiative model can offer a better description of the compiled data points that sample later epochs (≥55 days). The detection of a relativistic jet during the early expansion phase could be indirectly inferred through the modeling of the light curves (e.g., Mattila et al. 2018). This is as the afterglow emission from the expanding shock resulting from the jet—CNM interaction can be substantiated by emission from internal shocks in the jet (e.g., Alexander et al. 2020). In the absence of clear signatures of relativistic beaming based on low brightness temperatures, the emission from an off-axis jet may be significantly suppressed at large inclination angles, implying a domination by the afterglow emission. If the jet can survive the passage through a dense CNM, or with a favorable line of sight that renders it off-axis with respect to the observer, it may be detected and resolved such as in the case of Arp 299B AT1 at late epochs (Mattila et al. 2018). The significantly larger distance of AT2019dsg (in relation to Arp 299B AT1), and the decreasing flux density during the VLBI epochs however pose impediments to the detection.

The total energy in the outflow is minimal when the energy density of the relativistic electrons accelerated by the forward shock and that of the magnetic field are in equipartition (Scott & Readhead 1977; Chevalier 1998). The peak flux density is assumed to be associated with the transition of the synchrotron self-absorption frequency ν_a through the observation band. The evaluated peak flux density and associated time are then used to estimate the minimal size, magnetic field strength, energy output, and ambient medium number density; all of these parameters are critical for a nonrelativistic expansion (Chevalier 1998; Anderson et al. 2020). The application involves a sensitive dependence on microphysical parameters including the fraction of the total energy density in the electron kinetic energy _e and magnetic field _B . We use (_e , _B ) of (1/3, 1/3) and (0.1, 0.01) as typical representatives of TDEs and AT2019dsg, respectively (e.g., Ho et al. 2019; Anderson et al. 2020; Cendes et al. 2021a).

The frequencies that mark transitions in the evolving spectrum and their time evolution are as follows (see Appendix D): the self-absorption frequency ν_a ∝ t^{−2α(p+5)/(p+4)}, the synchrotron frequency ν_m ∝ t^−2α , and the synchrotron cooling frequency ν_c ∝ t^6α−2. The ordering ν_m ≤ ν_a as expected for a nonrelativistic outflow requires that the microphysical parameters _e = 0.1 and _B ≤ 0.01 at t = t_p . With these dependencies and using the above microphysical parameters, the frequencies at t = t_p are derived as ν_a = 4.59 × 10⁹ Hz, ν_m ≤ 2.56 × 10⁹ Hz and ν_c ≥ 1.92 × 10¹² Hz. This suggests that a spectral regime change to ν_m < ν_a occurred well before or near the peak of the light curve from an initial ordering ν_a < ν_m < ν_c with the observation frequency ν transitioning to the optically thin regime immediately after. The minimal values of the equipartition size are R_eq = 4.7 × 10¹⁶ cm (∼3 × 10⁴ R_S , where R_S = 1.59 × 10¹² cm is the Schwarzschild radius of a 5.4 × 10⁶ M_⊙ SMBH), magnetic field strength B_eq = 0.17 G, and energy output E_eq = 4.9 × 10⁴⁹ erg, and the ambient number density is n_eq = 5.7 × 10³ cm⁻³ (see Appendix D), indicating an energetic TDE (compared to thermal dominated TDEs, e.g., Alexander et al. 2020). As the outflow is inferred to accelerate in the early phase (Stein et al. 2021), it can have a constant velocity in a short phase between the flux-density peak (t_p = 152.8 days) and the first VLBI epoch (200 days), with a lower limit of R_eq/t_p ≈ 0.12 c for t_p = 152.8 days.

The partial disruption of a star (50%) with mass M_⋆ can potentially produce an energy E_⋆ = η(M_⋆/2)c². Assuming an accretion conversion efficiency η = 0.01 and equating with the energy output in the outflow (Equation (F10) with a SMBH mass M_• = 5.4 × 10⁶ M_⊙), we obtain a stellar mass in the range 0.14–2.05 M_⊙ for a penetration factor β in the range $1\mbox{--}7.86{M}_{\bullet }^{-2/3}{({M}_{\star }/{M}_{\odot })}^{-0.07}$ (see Appendix E). The corresponding β is 1–2.9, the outflow energy is ≤1.8 × 10⁵² erg, and the initial outflow velocity is 0.25–0.39 c, consistent with the expected deceleration to 0.12 c. The initial mass-outflow rate is (3.1–4.9) × 10⁻² M_⊙ yr⁻¹ (for t = t_p ). The slowdown is likely enabled by the outflow interaction with the surrounding CNM, owing to a density contrast ratio between the CNM and outflow of ≥32.3 at t = t_p indicating a denser CNM (using Equation (F12)). This could be aided by contributions from radiation drag, especially for high accretion rates and in the vicinity of the SMBH, which acts to make the outflow radially directed and at a saturated constant velocity (e.g., Mohan & Mangalam 2015; Yang et al. 2021).

The outflow can originate either from a wind driven by the radiation pressure of the accretion disk (e.g., Mohan & Mangalam 2015; Yang et al. 2021) or from unbound debris from the stellar disruption (e.g., Krolik et al. 2016). In the former scenario, the lifetime is between 40.9 and 228.0 days, set by the requirement that the shrinking photospheric radius equals the accretion-disk size (≈10 R_S , R_S is the Schwarzschild radius of the black hole), using Equation (F8). The corresponding disk luminosity is ∼0.2 times the Eddington luminosity, which can support nonrelativistic winds with velocities ≤0.2 c driven by radiation pressure from the disk emission (e.g., Yang et al. 2021). In the latter scenario, comparing the outflow velocity inferred here (≈0.12 c) to that associated with the outflow (Krolik et al. 2016) yields a nonphysical stellar mass of ≥1.3 × 10¹³ M_⋆. An origin as radiation-driven winds is thus more likely.

The acceleration of ultra-high-energy cosmic rays (e.g., Farrar & Gruzinov 2009) is closely associated with mechanisms of neutrino production (e.g., Guépin et al. 2018). This is naturally enabled in the case of a relativistic jet, through hadronic and photohadronic interactions (e.g., Biehl et al. 2018; Guépin et al. 2018). In the proposed concordance scenario (Winter & Lunardini 2021), the neutrino production is attributed to the interaction between X-ray photons (back-scattered from an expanding optically thick cocoon) and protons accelerated by internal shocks in a relativistic jet. As the density contrast is not especially large, a dark or hidden relativistic jet appears unlikely, also owing to the nondetection of a significant proper motion. For a relativistic jet with a bulk Lorentz factor equal to the Doppler boosting factor, the apparently superluminal velocity v_⊥ is related to the jet viewing angle θ as $({v}_{\perp }/c)=1/\tan \theta$. For a θ = 10°–45°, v_⊥ = 1–5.67 c corresponds to a proper motion of 0.29–1.67 mas yr⁻¹, which differs from our nondetection, thus rendering the proposed scenario for neutrino production by photohadronic interactions in a relativistic jet less likely.

We then explore nonjetted scenarios enabling the production of neutrinos. These include an association with the outflow or the accretion process (e.g., Hayasaki & Yamazaki 2019; Murase et al.2020), both of which are potential candidates, but neither can achieve the expected all-neutrino fluence (Murase et al. 2020). A detailed elaboration of the neutrino production is beyond the current scope. We thus discuss the relative plausibility of these scenarios broadly by estimating the maximum energy to which the protons (constituting the cosmic rays that produce the neutrinos through various reaction channels) can be accelerated and the associated timescales. In the former (outflow), the protons are accelerated by the outflow and the maximum energy (e.g., Zhang et al. 2017) is E_p ∝ t⁻¹ ≤ 1.3 × 10¹⁶ eV (at a time t = 175 days), limited by the dynamic timescale (see Equation (F13)). The detected neutrino energy E_ν ≈ 1.5% of E_p indicating that this is a possible channel. In the latter (accretion), the super-Eddington regime (with a lifetime t_j ≈ 35.3–576.3 days) that follows the debris circularization (≤11.1 days from Equation (F2)) is found to be inefficient in producing relativistic protons as the cooling from proton-proton collisions effectively limits the maximum energy (Hayasaki & Yamazaki 2019). In the radiatively inefficient accretion regime, the acceleration of protons by the plasma turbulence is limited by their diffusion with a consequent maximum energy E_p ≤ 3.7 × 10¹⁶ eV (see Equation (F14)) with the neutrino energy being ≈0.5% of E_p . However, the timescale for the transition from the super-Eddington to the radiatively inefficient accretion flow is long (≈16t_j ≥ 5.65 × 10² days), rendering it unlikely to occur during the current observation period. We thus tentatively favor the outflow scenario for the production of neutrinos.

The neutrino production itself may be enabled by either the p γ or the pp interactions at the base of the outflow, taken to be the photospheric radius, R_ph. By equating R_ph from Equation (F8) with an accretion-disk size ≈10 R_S (photosphere outside but in the vicinity of the disk), the resultant solution of R_ph/(10 R_S ) ≈ 1.56 constrains the stellar mass M_⋆ in the range 1.36–2.05 M_⊙ with a penetration factor β = 1 indicating that the disruption involves an extremely close approach to the SMBH. The photosphere is just outside and likely to be spatially distributed (e.g., Strubbe & Quataert 2009) around the disk. This site provides a suitable cross section for the interactions with a particle number density ≥10⁸ cm⁻³ (Stein et al. 2021) for γ ≥ 1.3 assuming that the CNM density scales as R^−γ in the innermost region. Furthermore, at the same epoch, the number density of the unbound stellar debris in the vicinity of the accretion disk is ≤3.9 × 10¹¹ cm⁻³ using Equation (F6) thus providing support to this scenario. A production at the radio emission region at R_eq ≈ 5.1 × 10¹⁶ cm is unlikely, owing to a relatively low particle number density n_eq = 4.1 × 10³ cm⁻³.

The VLBI observations cover the declining phase of the source evolution, and thus provide unique constraints on the flux density and proper motion. These play a key role in inferring a decelerating nonrelativistic expansion, and in reducing the free parameters in the radiative model. Similar to CNSS J0019+00 and ASASSN-14li, the relatively radio-quiet (≈4 × 10³⁸ erg s⁻¹) nature of AT2019dsg with a predominantly thermal emission and the nondetection of a relativistic jet are consistent with the nonthermal radio emission being from outflow activity. The large density contrast (≥32.3), albeit during the later epochs, hints at an outflow that may have decelerated owing to a frustration by the CNM (Stein et al.2021).

Using the results from the present study together with the measured flux densities during the rising phase in the literature, the following scenario emerges: a star of ≈1.35–2.05 M_⊙ is disrupted; a nonrelativistic outflow with a velocity of ≈0.1 c, originating from and radiatively driven by gas accretion onto the SMBH, interacts with the surrounding CNM at ≥4.7 × 10¹⁶ cm (∼3 × 10⁴ R_S ); the interactions result in an expanding forward shock that accelerates ambient electrons (constituting the CNM) to relativistic energies, thus producing the observed synchrotron radio emission; the production of neutrinos is likely at the base of the outflow, where protons can be efficiently accelerated to Peta eV energies and provide suitable cross sections for hadronic or photohadronic interactions, which is consistent with such an expectation from the detailed analysis in Stein et al. (2021). The origin of neutrinos inferred here is then in conflict with a possible origin from internal shocks in a relativistic jet (Winter & Lunardini 2021) or even a nonexpectation as espoused by Cendes et al. (2021a).

As the VLBI observations offer high resolutions, they offer a complementary (compared to multiwavelength electromagnetic observations) but unique perspective on the source evolution. Studies of TDEs, especially those nearby, are crucial to probe the outflow nature and constrain the CNM properties of the galactic nuclear regions. It is crucial that the VLBI monitoring of TDEs in general be initiated during the early evolutionary phases of the source evolution. This will help capture any possible relativistic expansion and the subsequent transition to a nonrelativistic phase. A sparse but active monitoring during later epochs can help identify the presence of a possible jet, if the conditions (nearby source, activity of the central engine, and density of the surrounding medium) are favorable for such an inference. This information is highly relevant for the accurate modeling of the afterglow radiative evolution that enables the extraction of physical parameters pertaining to the TDE physical and geometrical properties. The inference of a nonrelativistic outflow, and its potential role in the production of neutrinos provides an expanded inventory of extragalactic sources of neutrinos other than the Galactic PeV accelerators (HESS Collaboration et al. 2016; Ge et al. 2021). This includes blazar jets (e.g., Mannheim 1993; Giommi et al. 2020), the dominant nonblazar AGN populations (e.g., Murase & Waxman 2016; Hooper et al. 2019), and slow transients including gamma-ray bursts (e.g., Rachen & Mészáros 1998; Waxman & Bahcall 1999), making multimessenger investigations of similar sources essential and exciting (e.g., Guépin & Kotera 2017). This could be relevant in planning strategies to enable the effective use of global facilities for multimessenger astronomy.

We thank the anonymous referee for the constructive comments and suggestions, the inclusion of which have improved the content and presentation of our work. This work is supported by National Key R&D Programme of China (grant number 2018YFA0404603) and NSFC (1204130). The European VLBI Network (EVN) is a joint facility of independent European, African, Asian, and North American radio astronomy institutes. Scientific results from data presented in this publication are derived from the following EVN project code(s): EM140, RSM04. We thank Wen Chen at Kunming Station for making the first observation. e-VLBI research infrastructure in Europe is supported by the European Union, Seventh Framework Programme (FP7/2007-2013) under grant agreement number RI-261525 NEXPReS. e-MERLIN is a national facility operated by the University of Manchester at Jodrell Bank Observatory on behalf of STFC. The VLBI data processing made use of the computing resources of the China Square Kilometre Array (SKA) Regional center prototype, funded by the Ministry of Science and Technology of China and the Chinese Academy of Sciences.

The 5 GHz European VLBI network (EVN) observations were carried out on 2019 October 27 and November 12, and 2020 February 28. The observational setup is presented in Table 1. It includes the EVN project code, date, frequency, time duration, bandwidth, and participating telescopes (which included a maximum of 20 VLBI telescopes). The longest baselines in the East-West (Yebes, Spain—Tianma, China) and North-South (Badary, Russia—Hartebeesthoek, South Africa) directions span up to ≈9000 and 9800 km, respectively, contributing to the mas-scale high resolution. The large number of participating stations provided a good sampling in the (u, v) plane, optimally covering the source and a resulting high-sensitivity image with a root mean square noise down to ≈0.015 mJy beam⁻¹.

All the observations were conducted by utilizing the phase-referencing technique (Beasley & Conway 1995). During the observations, J2052+1619 served as the primary phase-referencing calibrator (Cal1, ∼200 mJy at 5 GHz), 2°.37 away from the target source AT2019dsg (Tar) on the plane of the sky. The nearby 20 mJy source J2058+1512 was used as the secondary phase calibrator (1°.08 away from the target). The bright (∼1 Jy) radio source J2031+1219 was observed as the fringe finder. The estimated VLBI flux densities of the above three calibrators are all obtained from the ASTROGEO Center, ⁵ where accurate astrometry information can also be found for selecting appropriate calibrators. The main observing cycle is "Cal1 (60s) - Tar (200s) - Cal1 (60s)." The secondary phase-referencing calibrator was observed for a scan of 200 s every five main cycles, serving as a reference to check the accuracy of the phase-referencing results. In sessions 1 and 3, the observational data at each station were recorded independently in disk modules and shipped to the correlator center hosted by the Joint Institute for VLBI ERIC (JIVE), located at Dwingeloo, the Netherlands. The total data recording rate is 2 gigabits per second. The disk packs from all stations were shipped to the EVN correlator SFXC (Keimpema et al. 2015). The correlation was carried out using the typical correlation setups for continuum sources (i.e., 1 s integration time and 1 MHz frequency resolution). In the second session (code: RSM04), the data recorded at each station were directly transferred in real time via broadband optical fiber to the SFXC software correlator in JIVE. Correlation was also carried out in the real-time mode. The rapid correlation and release of the second-session data enable us to quickly obtain the radio image of AT2019dsg 216 days post burst, validating the successful detection from the VLBI observation.

The correlated data were calibrated by following the standard procedures using the NRAO Astronomical Image Processing System (AIPS) software (Greisen 2003, p. 109). After loading the data with task FITLD, the interferometric visibility amplitudes were calibrated using the system temperatures and antenna gain curves measured at each telescope using task ANTAB. For some antennas (EM140A: Km; RSM04: Bd, Ir, Da, Kn, and Pi) without system temperature records during the observations, we initially used the nominal values of the system-equivalent flux density presented in the EVN status table ⁶ instead. The beginning 30 s of each scan were flagged due to the occupation of antenna slewing. The ionospheric dispersive delays were corrected using a map of the total electron content provided by Global Positioning System satellite observations ⁷ and applied to the visibilities through task TECOR. The VLBA procedure task VLBAPANG was used to correct phase variations caused by parallactic angles. A manual phase calibration and bandpass calibration were carried out using the best scans of the fringe finder through tasks FRING and BPASS. The global fringe-fitting was then performed on the phase-referencing calibrator J2052+1619.

After fringe-fitting of the primary phase-referencing calibrator, their visibility data were exported from AIPS and imported into Caltech Difmap software package (Shepherd 1997) for self-calibration and imaging. For these antennas, direct use of the nominal T_sys values could result in a large uncertainty on the source amplitude calibration. To avoid this issue, we first built CLEAN models of the phase-referencing calibrator only using antennas with T_sys measurements. After good gain solutions were obtained, we fixed the amplitude gain factors of these antenna by setting them to the reference antenna; next, we ran the amplitude and phase self-calibration again for all antennas (i.e., including those without T_sys values). Solving the closure amplitude then yielded the amplitude scaling factors for the antennas with nominal T_sys values. The corresponding antenna gain corrections from the amplitude self-calibration were then applied to the data using AIPS task CLCOR. Additional iteration of fringe-fitting was performed on the phase-referencing calibrator by taking account of its CLEAN component model to remove the small phase variations due to the core-jet structure of the phase calibrator. Then the obtained solutions were finally interpolated and applied to the secondary calibrator and the target source using CLCAL. The calibrated data were split into single-source files and imported into Difmap for imaging and model-fitting. The above data reduction and imaging procedures are the same for all the sessions. During the data processing, we found that some antennas have relatively large phase errors due to the coherent loss or bad weather; we thus excluded these antennas (EM140A: Sv, Zc, Hh, and Wb; RSM04: Zc, Hh, Wb, Cm, and De; EM140B: Sv, Zc, Hh, and Wb). The participating antennas are mentioned in Table 1. Given the weakness of the target source, no self-calibration was conducted on it. After CLEANing the peak regions of the map (i.e., signal-to-noise level ≥5σ), the natural grid weighting was used to create the final images. We used the Modelfit procedure in Difmap to constrain the observing parameters of the source (e.g., positions, sizes, and flux densities). Circular Gaussian models were used in the model fitting. Table 3 presents the fitting results. Uncertainties in the parameters are calculated based on those expected from an imaging analysis (Fomalont 1999) and include additional inputs. For the peak intensity and integrated flux density, an extra uncertainty of 5% was coupled with the original errors. The size errors are derived by using the ratio of the minor FWHM size of the synthesised beam to the component signal-to-noise ratio (S/N; peak intensities divided by rms noises). The following equations are used to estimate the above uncertainties (e.g., Fomalont 1999):

$$\begin{eqnarray}&&{\sigma }_{\mathrm{peak}}=\sqrt{{\sigma }_{\mathrm{rms}}^{2}+{(0.05\,{S}_{\mathrm{peak}})}^{2}},\end{eqnarray} \tag{ A1 }$$

$$\begin{eqnarray}&&{\sigma }_{\mathrm{tot}}=\sqrt{{\sigma }_{\mathrm{rms}}^{2}+{\left(\frac{{S}_{\mathrm{tot}}}{S/N}\right)}^{2}+{(0.05\,{S}_{\mathrm{tot}})}^{2}},\end{eqnarray} \tag{ A2 }$$

$$\begin{eqnarray}&&{\sigma }_{\theta }=\frac{{\theta }_{\mathrm{FWHM}}}{S/N},\end{eqnarray} \tag{ A3 }$$

where σ_rms is the rms noise of the images and S/N is the signal-to-noise ratio $\left(\displaystyle \frac{{S}_{\mathrm{peak}}}{{\sigma }_{\mathrm{rms}}}\right)$, and σ_peak, σ_tot, and σ_θ are the uncertainties of the peak flux density, the integrated flux density, and the FWHM size of a fitted Gaussian component, respectively.

The first full-array 8 hr EVN observation was conducted on 2019 October 27 (200 days after TDE). It was proposed as part of a two-epoch monitoring program of AT2019dsg to monitor the flux density and size evolution of the afterglow and to infer possible proper motion. A compact source with an integrated flux density of 0.61 mJy was constrained to a size ≤0.17 mas.

Following this, a target of opportunity of a 3 hr e-EVN observation was conducted on 2019 November 12 (216 days after the TDE). The aim was to verify the compactness of the source and the potential for long-term monitoring of the source. A compact source was successfully detected from this fast-response observation with an integrated flux density of 0.674 ± 0.048 mJy. The peak flux density is ≈28 times the rms noise, providing reliable detection confidence in the results from the first EVN session and detection promise for the next full-disk observation session scheduled three months later.

The second full-array 8 hr EVN observation was carried out on 2020 February 28 (session 3; 304 days after TDE), mainly to infer any possible proper motion, if the source remained detectable. A compact source was again detected, although its size expanded to ≤0.48 mas, and its integrated flux density was reduced to 0.44 mJy.

The astrometric position derived by averaging the peak positions over the three epochs is α(J2000) =20^h57^m0296434 ± 000001 and δ(J2000) = 14°1'16''.29425 ±0''.00020. The uncertainties are the square root of the sum of the squares of two error terms: the astrometric error of the phase-referencing calibrator Cal1 (0.13 mas on R.A. and 0.17 mas on decl., from ASTROGEO) and the simulated astrometric accuracy of VLBI phase-referenced observations (that gives a statistical positional error of Δ_R.A. ≈ 0.06 mas and Δ_decl. ≈ 0.10 mas; see Pradel et al. 2006). The total uncertainty is dominated by the astrometric accuracy of the calibrator.

Errors on the measured source positions are composed of systematic and random contributions. Systematic errors include contributions from instrumental and imaging sources. We followed a similar approach to that of Mohan et al. (2020), with the employment of two phase calibrators CAL1 (J2052+1619) and CAL2 (J2058+1512). To minimize systematic errors, in each epoch, we used the calibration solutions from CAL1 to derive the positions for both CAL2 and the target source (AT2019dsg). Then we calculated the relative positions between AT2019dsg and CAL2 (see column 5 in Table 2); this helped monitor positional changes across the three observation epochs. Since CAL2 is an extragalactic radio quasar, and the positions used correspond the the bright radio core, CAL2 can be treated as a stable source across the three epochs. The above methodology renders most instrumental- and calibration-based errors negligible. The remaining contributors to systematic errors are then only from the measurement errors on CAL2 and AT2019dsg, as well as the phase-referencing errors (i.e., the error caused by the separation between AT2019dsg and the phase calibrator; see Pradel et al. 2006 for details). The measurements and errors are listed in Table 2. As we use the relative position (w.r.t. CAL2) to identify any shift in the emission center of AT2019dsg, a major portion of the systematic errors can be eliminated. The systematic errors σ_sys reported in column 6 of Table 2 are estimated as the standard deviation on the CAL2 offsets during the three observation epochs. For the random error σ_ran, a fitting error similar to that suggested in Fomalont (1999) was applied, with σ_ran = θ_FWHM/(2 S/N). For resolved components, the θ_FWHM are the fitted component sizes; for unresolved sources, θ_FWHM are the synthesised beams. In columns 7, the random errors for both the target and CAL2 were involved. The final errors (columns 8) were calculated by combining the error from both contributions (columns 6 and 7). Using these measurements (columns 5 and 8 of Table 2), the relative position change of AT2019dsg w.r.t. CAL2 is found to lie at a midpoint (RA, DEC) = (0.06 ± 0.12 mas, −0.07 ± 0.13 mas) (1σ uncertainties; see Figure 1 panel 4, which illustrates these measurements and the 3σ uncertainty ellipse enclosing them) indicating that there is no significant shift in the emission center of AT2019dsg across the three observational epochs. The distance between the midpoint and any of the other two points sets a maximum relative position change of 0.32 ± 0.22 mas (1σ).

The compiled 5 GHz light curve (see Figure 2) includes data points covering the rising phase and peak (Stein et al. 2021; Cannizzaro et al. 2021), and our three VLBI points in the declining phase. Data points from Cendes et al. (2021a) are not employed in the fit owing to their VLA flux densities being systematically lower than the e-MERLIN measurements reported in Cannizzaro et al. (2021). As our VLBI estimates are likely to be closer to the e-MERLIN values, we mainly use the former measurements only when conducting the light-curve fitting.

The light curve is fitted using the weighted nonlinear least squares method with a model of the form ${F}_{\nu }={F}_{{\nu }_{0}}+{F}_{{\nu }_{1}}{({(t/{t}_{p})}^{5{\alpha }_{1}}+{(t/{t}_{p})}^{5{\alpha }_{2}})}^{-1/5}$, with weights based on the flux-density measurement errors, the flux-density normalization $0.0\lt {F}_{{\nu }_{0}}\leqslant 0.5$ mJy, $0.0\lt {F}_{{\nu }_{1}}\leqslant 3.0$ mJy, the rising slope 0.0 ≤ α₁ ≤ 5.0 and declining slope − 10.0 ≤ α₂ ≤ 0.0, a turnover time t_p days, and a smoothness parameter fixed at 5. The availability of flux-density measurements during the last three epochs based on our EVN observations plays a crucial role in constraining the optically thin declining phase. The fit yields the parameters ${F}_{{\nu }_{0}}=0.17\pm 0.17$ mJy, ${F}_{{\nu }_{1}}=1.18\pm 0.18$ mJy, α₁ = 2.38 ± 0.93, and α₂ = − 2.56 ± 1.31. These yield a peak flux density ${F}_{{\nu }_{p}}=1.19\pm 0.18$ mJy and an associated peak time t_p = 152.80 ± 16.19 days. The weighted nature of the fitting through the dependence on measurement errors in the data points is reflected in the relatively large uncertainties above. This is especially the case near the light-curve peak (flux density and time of occurrence) and in the declining slope where a paucity of data points and measurement errors render larger imprints. The fit parameters are summarized in Table 4.

The interaction of the outflow with the CNM produces an expanding forward shock. Assuming an initial relativistic expansion, the Lorentz factor of the shock Γ_sh and the projected size r of the emitting region are often assumed to follow a self-similar evolution (Blandford & McKee 1976), and are given by (e.g., Metzger et al. 2012)

$$\begin{eqnarray}&&{{\rm{\Gamma }}}_{\mathrm{sh}}={\left(\displaystyle \frac{17-4k}{16\pi }\displaystyle \frac{{E}_{j}}{{m}_{p}{{nr}}^{3}{c}^{2}}\right)}^{1/2},\end{eqnarray} \tag{ C1 }$$

$$\begin{eqnarray}&&r=\displaystyle \frac{2{{\rm{\Gamma }}}_{\mathrm{sh}}^{2}{ct}}{(1+z)},\end{eqnarray} \tag{ C2 }$$

where the number density is assumed to scale with the distance as $n={n}_{0}{(r/{r}_{0})}^{-k}$ for a scaling constant k = 0 in the case of a surrounding medium of a constant number density, and k = 2 for a stratified medium, for a fiducial distance r₀ with associated number density n₀. The total kinetic energy of the outflow is E_j . Using the scalings E_j,52 = E_j /(10⁵² erg s⁻¹), t_d = t/(86400 s), and D_L,100 = D_L /(100 Mpc) and assuming r₀ = 1 pc, the parameters Γ_sh and θ_A = r/D_A (where θ_A is the source angular size and D_A = D_L /(1 + z)² is the source angular distance) are

$$\begin{eqnarray}&&{{\rm{\Gamma }}}_{\mathrm{sh}}{\left(k=0)=7.96\,{\left(\displaystyle \frac{{E}_{j,52}}{{n}_{0}}\right)}^{1/8}{t}_{d}^{-3/8}(1+z\right)}^{3/8},\end{eqnarray} \tag{ C3 }$$

$$\begin{eqnarray}&&{\theta }_{A}{\left(k=0)=(0.18\,\mathrm{mas})\,{\left(\displaystyle \frac{{E}_{j,52}}{{n}_{0}}\right)}^{1/4}{t}_{d}^{1/4}(1+z\right)}^{7/4}{D}_{L,100}^{-1},\end{eqnarray} \tag{ C4 }$$

$$\begin{eqnarray}&&{{\rm{\Gamma }}}_{\mathrm{sh}}{\left(k=2)=2.21\,{\left(\displaystyle \frac{{E}_{j,52}}{{n}_{0}}\right)}^{1/4}{t}_{d}^{-1/4}(1+z\right)}^{1/4},\end{eqnarray} \tag{ C5 }$$

$$\begin{eqnarray}&&{\theta }_{A}{\left(k=2)=(0.02{\rm{mas}}){\left(\displaystyle \frac{{E}_{j,52}}{{n}_{0}}\right)}^{1/2}{t}_{d}^{1/2}(1+z\right)}^{3/2}{D}_{L,100}^{-1}.\end{eqnarray} \tag{ C6 }$$

We assume that the radio afterglow emission originates from a nonrelativistic, spherically expanding forward shock. In this formulation (Chevalier 1998; Anderson et al. 2020), the peak of the light curve corresponds to the transition of the synchrotron emission from optically thick to thin through self-absorption. The relativistic electrons responsible for the radio emission are assumed to have a power-law distribution in energy E, with the number per unit of energy N = N₀ E^−p where N₀ is a constant evaluated at a given time and p ≥ 2 is the power-law index.

The flux densities in the optically thick ( ∝ ν^5/2) and thin ( ∝ ν^−(p−1)/2) regimes are

$$\begin{eqnarray}&&{F}_{\nu }=\displaystyle \frac{\pi {R}^{2}}{{D}_{L}^{2}}\displaystyle \frac{{c}_{5}}{{c}_{6}}{B}^{-1/2}{\left(\displaystyle \frac{\nu }{2{c}_{1}}\right)}^{5/2},\end{eqnarray} \tag{ D1 }$$

$$\begin{eqnarray}&&{F}_{\nu }=\displaystyle \frac{4\pi {R}^{3}f}{3{D}_{L}^{2}}{c}_{5}{N}_{0}{B}^{(p+1)/2}{\left(\displaystyle \frac{\nu }{2{c}_{1}}\right)}^{-(p-1)/2},\end{eqnarray} \tag{ D2 }$$

where R is the size of the emitting region, D_L is the luminosity distance, B is the magnetic field strength, and f is the emission filling factor. The radiative constants c₁, c₅, and c₆ as functions of p are given by Pacholczyk (1970):

$$\begin{eqnarray}&&{c}_{1}=\displaystyle \frac{3e}{4\pi {m}_{e}^{3}{c}^{5}},\end{eqnarray} \tag{ D3 }$$

$$\begin{eqnarray}&&{c}_{5}=\displaystyle \frac{{3}^{1/2}{e}^{3}}{4\pi {m}_{e}{c}^{2}(p+1)}{\rm{\Gamma }}\left(\displaystyle \frac{p}{4}+\displaystyle \frac{19}{12}\right){\rm{\Gamma }}\left(\displaystyle \frac{p}{4}-\displaystyle \frac{1}{12}\right),\end{eqnarray} \tag{ D4 }$$

$$\begin{eqnarray}&&{c}_{6}=\displaystyle \frac{{3}^{1/2}{e}^{3}}{8\pi {m}_{e}}{\left(\displaystyle \frac{3e}{2\pi {m}_{e}^{3}{c}^{5}}\right)}^{-2}{\rm{\Gamma }}\left(\displaystyle \frac{p}{4}+\displaystyle \frac{1}{6}\right){\rm{\Gamma }}\left(\displaystyle \frac{p}{4}+\displaystyle \frac{11}{6}\right)\end{eqnarray} \tag{ D5 }$$

where the electric charge e = 4.8 × 10⁻¹⁰ e.s.u., electron mass m_e = 9.11 × 10⁻²⁸ g, speed of light c = 3 × 10¹⁰ cm s⁻¹, and Γ represents the gamma function. The constant energy density term N₀ in Equation (D2) is obtained from the equipartition condition

$$\begin{eqnarray}&&{N}_{0}=\left(\displaystyle \frac{{\epsilon }_{e}}{{\epsilon }_{B}}\right)\displaystyle \frac{{B}^{2}}{8\pi }(p-2){E}_{l}^{p-2},\end{eqnarray} \tag{ D6 }$$

where _e and _B are the fractions of the total energy density in the particle kinetic energy density and in the magnetic field respectively, and E_l = 0.51 MeV = 8.2 × 10⁻⁷ erg is the electron rest mass-energy and is representative of the energy at which the electrons transition to relativistic.

Equating the flux densities from Equations (D1) and (D2), and using Equation (D6), the minimal size R_eq and magnetic field strength B_eq are expressed in terms of the peak flux density ${F}_{{\nu }_{p}}$ and corresponding frequency ν_p as

$$\begin{eqnarray}\begin{array}{rcl}{R}_{\mathrm{eq}} & = & {\left(\displaystyle \frac{6{c}_{6}^{p+5}{F}_{{\nu }_{p}}^{p+6}{D}_{L}^{2p+12}}{({\epsilon }_{e}/{\epsilon }_{B})f(p-2){\pi }^{p+5}{c}_{5}^{p+6}{E}_{l}^{p-2}}\right)}^{1/(2p+13)}{\left(\displaystyle \frac{{\nu }_{p}}{2{c}_{1}}\right)}^{-1}\\ & = & (1.13\times {10}^{17}\,\mathrm{cm}){\left(\displaystyle \frac{{F}_{{\nu }_{p}}}{\mathrm{mJy}}\right)}^{0.47}\\ & & \times {\left(\displaystyle \frac{{D}_{L}}{100\,\mathrm{Mpc}}\right)}^{0.95}{\left(\displaystyle \frac{{\epsilon }_{e}}{{\epsilon }_{B}}\right)}^{-0.05}{f}^{-0.05}{\left(\displaystyle \frac{{\nu }_{p}}{\mathrm{GHz}}\right)}^{-1},\end{array}\end{eqnarray} \tag{ D7 }$$

$$\begin{eqnarray}\begin{array}{rcl}{B}_{\mathrm{eq}} & = & {\left(\displaystyle \frac{36{\pi }^{3}{c}_{5}}{{({\epsilon }_{e}/{\epsilon }_{B})}^{2}{f}^{2}{(p-2)}^{2}{c}_{5}^{3}{E}_{l}^{2(p-2)}{F}_{{\nu }_{p}}{D}_{L}^{2}}\right)}^{2/(2p+13)}\left(\displaystyle \frac{{\nu }_{p}}{2{c}_{1}}\right)\\ & = & (0.06\,{\rm{G}}){\left(\displaystyle \frac{{F}_{{\nu }_{p}}}{\mathrm{mJy}}\right)}^{-0.11}{\left(\displaystyle \frac{{D}_{L}}{100\,\mathrm{Mpc}}\right)}^{-0.22}\\ & & \times {\left(\displaystyle \frac{{\epsilon }_{e}}{{\epsilon }_{B}}\right)}^{-0.22}{f}^{-0.22}\left(\displaystyle \frac{{\nu }_{p}}{\mathrm{GHz}}\right),\end{array}\end{eqnarray} \tag{ D8 }$$

and the minimal energy and number density are evaluated as

$$\begin{eqnarray}\begin{array}{rcl}{E}_{\mathrm{eq}} & = & \displaystyle \frac{{B}_{\mathrm{eq}}^{2}}{8\pi {\epsilon }_{B}}\displaystyle \frac{4\pi }{3}{R}_{\mathrm{eq}}^{3}f\\ & = & (1.11\times {10}^{48}\,\mathrm{erg}){\left(\displaystyle \frac{{F}_{{\nu }_{p}}}{\mathrm{mJy}}\right)}^{1.20}{\left(\displaystyle \frac{{D}_{L}}{100\,\mathrm{Mpc}}\right)}^{2.41}\\ & & \times {\left(\displaystyle \frac{{\epsilon }_{e}}{{\epsilon }_{B}}\right)}^{-0.59}{\epsilon }_{B}^{-1}{f}^{0.41}{\left(\displaystyle \frac{{\nu }_{p}}{\mathrm{GHz}}\right)}^{-1},\end{array}\end{eqnarray} \tag{ D9 }$$

$$\begin{eqnarray}\begin{array}{rcl}{n}_{\mathrm{eq}} & = & \displaystyle \frac{{\epsilon }_{e}}{{\epsilon }_{B}}\displaystyle \frac{{B}_{\mathrm{eq}}^{2}}{8\pi }\left(\displaystyle \frac{p-2}{p-1}\right){E}_{l}^{-1}\\ & = & (91.76\,{\mathrm{cm}}^{-3}){\left(\displaystyle \frac{{F}_{{\nu }_{p}}}{\mathrm{mJy}}\right)}^{-0.22}{\left(\displaystyle \frac{{D}_{L}}{100\,\mathrm{Mpc}}\right)}^{-0.44}\\ & & \times {\left(\displaystyle \frac{{\epsilon }_{e}}{{\epsilon }_{B}}\right)}^{0.56}{f}^{-0.44}{\left(\displaystyle \frac{{\nu }_{p}}{\mathrm{GHz}}\right)}^{2},\end{array}\end{eqnarray} \tag{ D10 }$$

where the expressions for R_eq, B_eq, E_eq, and n_eq are evaluated for p = 2.7 and assuming fiducial values for the free parameters.

The analysis assumes that the peak frequency ν_p is equal to the synchrotron self-absorption frequency ν_a , which is obtained from the condition that the optical depth to self-absorption τ ≈ α_ν R = 1 at ν = ν_a , where α_ν is the absorption coefficient, and R is the time-evolving size of the emission region. The time-dependent physical parameters include $R\propto {t}^{\alpha }={R}_{\mathrm{eq}}{(t/{t}_{p})}^{\alpha }$, $B\propto {R}^{-2}\propto {t}^{-2\alpha }={B}_{\mathrm{eq}}{(t/{t}_{p})}^{-2\alpha }$ and ${N}_{0}\propto {B}^{2}={N}_{0,{eq}}{(t/{t}_{p})}^{-4\alpha }$ where N_0,eq = N₀(B = B_eq) (see Equation (D6)); t_p is the time associated with the light-curve peak.

The absorption coefficient is

$$\begin{eqnarray}&&{\alpha }_{\nu }={c}_{6}{N}_{0}{B}^{(p+2)/2}{\left(\displaystyle \frac{\nu }{2{c}_{1}}\right)}^{-(p+4)/2}.\end{eqnarray} \tag{ D11 }$$

Using the condition α_ν R = 1 and the time-dependent parameters, the synchrotron self-absorption frequency is

$$\begin{eqnarray}\begin{array}{rcl}{\nu }_{a} & = & 2{c}_{1}{({{Rc}}_{6}{N}_{0})}^{2/(p+4)}{B}^{(p+2)/(p+4)}\propto {t}^{-2\alpha (p+5)/(p+4)}\\ & = & (0.92\times {10}^{9}\,\mathrm{Hz}){f}^{-0.30}\left(\displaystyle \frac{{\nu }_{p}}{\mathrm{GHz}}\right){\left(\displaystyle \frac{t}{{t}_{p}}\right)}^{-2.30\alpha }.\end{array}\end{eqnarray} \tag{ D12 }$$

The characteristic synchrotron frequency emitted by a relativistic electron with a Lorentz factor γ_e is

$$\begin{eqnarray}&&{\nu }_{e}=\displaystyle \frac{{eB}{\gamma }_{e}^{2}}{2\pi {m}_{e}c}.\end{eqnarray} \tag{ D13 }$$

For ${\gamma }_{e}={\gamma }_{m}=\left(\tfrac{p-2}{p-1}\right)\left(\tfrac{{m}_{p}}{{m}_{e}}\right){\epsilon }_{e}$, the minimum injected Lorentz factor, the corresponding synchrotron frequency is

$$\begin{eqnarray}&&\begin{array}{l}{\nu }_{m}\,=\,\frac{{eB}{\gamma }_{m}^{2}}{2\pi {m}_{e}c}\propto {t}^{-2\alpha }\\ \,=\,(1.08\times {10}^{11}\,\mathrm{Hz}){\left(\frac{{F}_{{\nu }_{p}}}{\mathrm{mJy}}\right)}^{-0.11}{\left(\frac{{D}_{L}}{100\,\mathrm{Mpc}}\right)}^{-0.22}\\ \times {\left(\frac{{\epsilon }_{e}}{{\epsilon }_{B}}\right)}^{-0.22}{\epsilon }_{e}^{2}{f}^{-0.22}\left(\frac{{\nu }_{p}}{\mathrm{GHz}}\right){\left(\frac{t}{{t}_{p}}\right)}^{-2\alpha };\end{array}\end{eqnarray} \tag{ D14 }$$

while for ${\gamma }_{e}={\gamma }_{c}=\tfrac{6\pi {m}_{e}c}{{\sigma }_{T}{B}^{2}t}$, the critical Lorentz factor above which synchrotron cooling dominates and where σ_T = 6.65 ×10⁻²⁵ cm² is the Thomson cross section, the corresponding synchrotron cooling frequency is

$$\begin{eqnarray}\begin{array}{rcl}{\nu }_{c} & = & \displaystyle \frac{18\pi {m}_{e}{ce}}{{\sigma }_{T}^{2}{B}^{3}{t}^{2}}\propto {t}^{6\alpha -2}\\ & = & (7.27\times {10}^{17}\,\mathrm{Hz}){\left(\displaystyle \frac{{F}_{{\nu }_{p}}}{\mathrm{mJy}}\right)}^{0.33}{\left(\displaystyle \frac{{D}_{L}}{100\,\mathrm{Mpc}}\right)}^{0.66}\\ & & \times {\left(\displaystyle \frac{{\epsilon }_{e}}{{\epsilon }_{B}}\right)}^{0.66}{f}^{0.66}{\left(\displaystyle \frac{{\nu }_{p}}{\mathrm{GHz}}\right)}^{-3}{\left(\displaystyle \frac{t}{{t}_{p}}\right)}^{6\alpha -2}{\left(\displaystyle \frac{{t}_{p}}{\mathrm{day}}\right)}^{-2}.\end{array}\end{eqnarray} \tag{ D15 }$$

Post the TDE, most gravitationally bound stellar material returns to the pericenter when orbiting the SMBH. The fallback timescale is (Lodato & Rossi 2011)

$$\begin{eqnarray}&&{t}_{\mathrm{fb}}=(41\,\mathrm{days})\,{M}_{6}^{1/2}{m}_{\star }^{(1-3\xi )/2}{\beta }^{-3},\end{eqnarray} \tag{ E1 }$$

assuming that the main-sequence mass–radius relationship ${R}_{\star }={R}_{\odot }{m}_{\star }^{1-\xi }$ (Kippenhahn & Weigert 1994) is valid for the disrupted star with a scaled mass m_⋆ = M_⋆/M_⊙ and radius R_⋆, where the index ξ ≈ 0.2 for 0.1 ≤ m_⋆ ≤ 1.0, and ξ ≈ 0.4 for m_⋆ ≥ 1.0 (Krolik & Piran 2012). The other scaled quantities include M₆ = M_•/(10⁶ M_⊙) and a penetration factor β ≡ R_t /R_p defined as the ratio of the tidal radius of the SMBH R_t to the pericenter distance R_p .

The tidal radius of the SMBH R_t can be expressed in terms of the SMBH mass and stellar properties (Evans & Kochanek 1989) as

$$\begin{eqnarray}&&{R}_{t}\approx {\left(\displaystyle \frac{{M}_{\bullet }}{{M}_{\star }}\right)}^{1/3}{R}_{\star }=(6.95\times {10}^{12}\,\mathrm{cm})\,{M}_{6}^{1/3}{m}_{\star }^{(2/3)-\xi },\end{eqnarray} \tag{ E2 }$$

which is then scaled in terms of the Schwarzschild radius R_S = 2GM_•/c² = (2.95 × 10¹¹ cm)M₆, as

$$\begin{eqnarray}&&({R}_{t}/{R}_{S})=(23.58)\,{M}_{6}^{-2/3}{m}_{\star }^{(2/3)-\xi }.\end{eqnarray} \tag{ E3 }$$

The distance of the closest approach to the SMBH is the pericenter distance R_p ≥ 3R_S , the innermost stable circular orbit for a nonrotating Schwarzschild black hole. With R_t ≥ R_p , R_t /R_S ≥ R_p /R_S ≥ 3, and the penetration factor β ≡ R_t /R_p is in the range $1\leqslant \beta \leqslant (7.86)\,{M}_{6}^{-2/3}{m}_{\star }^{(2/3)-\xi }$.

The fallback accretion rate onto the SMBH (Evans & Kochanek 1989) is

$$\begin{eqnarray}&&{\dot{M}}_{\mathrm{fb}}=\displaystyle \frac{1}{3}\displaystyle \frac{{M}_{\star }}{{t}_{\mathrm{fb}}}{\left(\displaystyle \frac{t}{{t}_{\mathrm{fb}}}\right)}^{-5/3},\end{eqnarray} \tag{ E4 }$$

at a time t post the tidal disruption. The peak rate corresponds to t = t_fb, where

$$\begin{eqnarray}&&{\dot{M}}_{{\rm{p}}}={\dot{M}}_{\mathrm{fb}}{| }_{t={t}_{\mathrm{fb}}}\approx (3.0\,{M}_{\odot }\,{\mathrm{yr}}^{-1})\,{M}_{6}^{-1/2}{m}_{\star }^{(1+3\xi )/2}{\beta }^{3}.\end{eqnarray} \tag{ E5 }$$

The Eddington accretion rate is

$$\begin{eqnarray}&&{\dot{M}}_{\mathrm{Edd}}=\displaystyle \frac{{L}_{\mathrm{Edd}}}{\epsilon {c}^{2}}=\displaystyle \frac{4\pi {{GM}}_{\bullet }{m}_{p}}{\epsilon c{\sigma }_{T}}=(0.023\,{M}_{\odot }{\mathrm{yr}}^{-1})\,{M}_{6},\end{eqnarray} \tag{ E6 }$$

for an assumed efficiency factor that we fix at a fiducial value = 0.1 and where σ_T = 6.65 × 10⁻²⁵ cm² is the Thomson cross section for electron scattering. Using the scaling t_d = t/(1 day), the fallback accretion rate and peak rate are super-Eddington (see Equations (E4) and Equation (E5)), with

$$\begin{eqnarray}&&{\dot{m}}_{\mathrm{fb}}=({\dot{M}}_{\mathrm{fb}}/{\dot{M}}_{\mathrm{Edd}})=(6.34\times {10}^{4})\,{M}_{6}^{-2/3}{m}_{\star }^{(4-3\xi )/3}{\beta }^{-2}{t}_{d}^{-5/3},\end{eqnarray} \tag{ E7 }$$

$$\begin{eqnarray}&&{\dot{m}}_{{\rm{p}}}=({\dot{M}}_{{\rm{p}}}/{\dot{M}}_{\mathrm{Edd}})=(130.43)\,{M}_{6}^{-3/2}{m}_{\star }^{(1+3\xi )/2}{\beta }^{3}.\end{eqnarray} \tag{ E8 }$$

During the super-Eddington accretion phase, the returning tidal stream post the stellar disruption can encounter the stellar debris being accreted onto the SMBH, resulting in a shocked debris and the subsequent launching of outflows (e.g., Jiang et al. 2016; Lu & Bonnerot 2020; Murase et al. 2020). The outflow is expected to be launched close to the circularization radius of the accreted material R_C ≈ 2R_p = 2R_t β⁻¹ with an initial velocity (e.g., Strubbe & Quataert 2009),

$$\begin{eqnarray}&&{v}_{w}\approx {\left(\displaystyle \frac{{{GM}}_{\bullet }}{{R}_{C}}\right)}^{1/2}=(0.1\,c)\,{M}_{6}^{1/3}{m}_{\star }^{-(1/3-\xi /2)}{\beta }^{1/2},\end{eqnarray} \tag{ F1 }$$

with the equality corresponding to the escape velocity at the radius R_C . Assuming that the efficiency for the mass-energy conversion during circularization is 1%, the associated timescale is (Hayasaki & Yamazaki 2019)

$$\begin{eqnarray}&&{t}_{\mathrm{circ}}=(1.8\times {10}^{6}\,{\rm{s}}){\beta }^{-3/2}{M}_{6}^{-1/2}{m}^{3(1-2\xi )/2}.\end{eqnarray} \tag{ F2 }$$

The fallback accretion rate declines as t^−5/3 (see Equation (E7)) with ${\dot{M}}_{\mathrm{fb}}={\dot{M}}_{\mathrm{Edd}}$ denoting a transition from the super-Eddington to the sub-Eddington regime. The corresponding timescale t_j is taken to be the lifetime of the outflow and is given by

$$\begin{eqnarray}&&{t}_{j}=(6.54\times {10}^{7}\,s)\,{M}_{6}^{-2/5}{m}_{\star }^{(4-3\xi )/5}{\beta }^{-6/5}.\end{eqnarray} \tag{ F3 }$$

With the above outflow velocity v_w and $n\propto {(R/{R}_{\mathrm{eq}})}^{-\gamma }={n}_{\mathrm{eq}}{(t/{t}_{p})}^{-\gamma \alpha }$ to represent the number density decline with distance $R={R}_{\mathrm{eq}}{(t/{t}_{p})}^{\alpha }$, the mass-outflow rate is

$$\begin{eqnarray}&&{\dot{M}}_{w}=4\pi {R}^{2}{v}_{w}{{nm}}_{p},\end{eqnarray} \tag{ F4 }$$

$$\begin{eqnarray}&&\begin{array}{l}=\ (1.18\times {10}^{-3}\,{M}_{\odot }\,{\mathrm{yr}}^{-1})\,{M}_{6}^{1/3}{m}_{\star }^{-(1/3-\xi /2)}{\beta }^{1/2}\\ \quad \times \,{\left(\displaystyle \frac{{F}_{{\nu }_{p}}}{\mathrm{mJy}}\right)}^{0.73}{\left(\displaystyle \frac{{D}_{L}}{100\,\mathrm{Mpc}}\right)}^{1.46}{\left(\displaystyle \frac{{\epsilon }_{e}}{{\epsilon }_{B}}\right)}^{0.46}{f}^{-0.54}{\left(\displaystyle \frac{t}{{t}_{p}}\right)}^{\alpha (2-\gamma )}.\end{array}\end{eqnarray} \tag{ F5 }$$

An upper limit on the particle number density in the vicinity of the accretion disk is based on the unbound stellar debris that is spatially and kinematically distributed in this region (Strubbe & Quataert 2009):

$$\begin{eqnarray}&&{n}_{a}\leqslant (1.46\times {10}^{18}\,{\mathrm{cm}}^{-3}){M}_{6}^{1/6}{m}_{\star }^{(5/6)-(3\xi /2)}{t}_{d}^{-3}.\end{eqnarray} \tag{ F6 }$$

Assuming that the constituent gas in the outflow is in thermal equilibrium with a consequent blackbody spectrum, the total luminosity is approximated by the Stefan–Boltzmann law (Strubbe & Quataert 2009) such that

$$\begin{eqnarray}\begin{array}{rcl}{L}_{j} & = & {\sigma }_{\mathrm{SB}}(4\pi {R}_{\mathrm{ph}}^{2}){T}_{\mathrm{ph}}^{4}=(6.68\times {10}^{44}\,\mathrm{erg}\,{{\rm{s}}}^{-1})\\ & & \times \,{M}_{6}^{8/9}{m}_{\star }^{(2-\xi )/6}{\beta }^{-2/3}{t}_{d}^{-5/9},\end{array}\end{eqnarray} \tag{ F7 }$$

where σ_SB = 5.67 × 10⁻⁵ erg cm⁻² s⁻¹ K⁻⁴ is the Stefan–Boltzmann constant; assuming a photospheric radius R_ph and an associated thermal temperature T_ph:

$$\begin{eqnarray}&&{R}_{\mathrm{ph}}=(6.3\times {10}^{16}\,\mathrm{cm})\,{m}_{\star }^{5/3-3\xi /2}{\beta }^{-2}{t}_{d}^{-5/3},\end{eqnarray} \tag{ F8 }$$

$$\begin{eqnarray}&&{T}_{\mathrm{ph}}=(1.25\times {10}^{3}\,{\rm{K}})\,{M}_{6}^{2/9}{m}_{\star }^{3/4-17\xi /24}{\beta }^{5/6}{t}_{d}^{25/36},\end{eqnarray} \tag{ F9 }$$

assuming a mass-outflow rate ${\dot{M}}_{w}\leqslant 0.1{\dot{M}}_{\mathrm{fb}}$ (e.g., Mageshwaran & Mangalam 2015) and the equality in Equation (F1). Also employed in the above equation are the expressions for R_p /(3R_S ) and ${\dot{m}}_{\mathrm{fb}}$ from Equation (E2) and Equation (E7), respectively. The corresponding total kinetic energy E_j is estimated using the above thermal luminosity and an additive component from the power-law tail during the declining phase:

$$\begin{eqnarray}\begin{array}{rcl}{E}_{j} & = & {\displaystyle \int }_{0}^{{t}_{j}}{L}_{j}{dt}+{\displaystyle \int }_{{t}_{j}}^{t}{L}_{j}{(t/{t}_{j})}^{-5/3}{dt}\\ & \approx & (3.81\times {10}^{51}\,\mathrm{erg})\,{M}_{6}^{32/45}{m}_{\star }^{(62-39\xi )/90}{\beta }^{-6/5}.\end{array}\end{eqnarray} \tag{ F10 }$$

The number density of material in the outflow is given by

$$\begin{eqnarray}&&{n}_{j}=\displaystyle \frac{{L}_{j}}{4\pi {R}^{2}{m}_{p}{c}^{3}}.\end{eqnarray} \tag{ F11 }$$

Using L_j from Equation (F7) and the scalings $n={n}_{\mathrm{eq}}{(t/{t}_{p})}^{-\gamma \alpha }$ and $R={R}_{\mathrm{eq}}{(t/{t}_{p})}^{\alpha }$, the density contrast between the CNM and the outflow,

$$\begin{eqnarray}\begin{array}{rcl}\kappa & = & \displaystyle \frac{n}{{n}_{j}}\approx (1.0)\,{M}_{6}^{-8/9}{m}_{\star }^{-(2-\xi )/6}{\beta }^{2/3}{\left(\displaystyle \frac{{F}_{{\nu }_{p}}}{\mathrm{mJy}}\right)}^{0.73}\\ & & \times {\left(\displaystyle \frac{{D}_{L}}{100\,\mathrm{Mpc}}\right)}^{1.46}{\left(\displaystyle \frac{{\epsilon }_{e}}{{\epsilon }_{B}}\right)}^{0.46}{f}^{-0.54}\\ & & \times {\left(\displaystyle \frac{t}{{t}_{p}}\right)}^{(2-\gamma )\alpha +5/9}{\left(\displaystyle \frac{{t}_{p}}{\mathrm{day}}\right)}^{5/9}.\end{array}\end{eqnarray} \tag{ F12 }$$

The maximum energy to which the protons can be accelerated by the outflow is limited by comparing the associated timescale with the dynamical timescale (e.g., Zhang et al. 2017). Using the scalings for the size, velocity, and magnetic field strength in the emission region, the maximum energy is

$$\begin{eqnarray}\,\begin{array}{rcl}{E}_{p} & \leqslant & \frac{3{v}^{2}}{20{c}^{2}}{eBct}=\frac{3}{20}\frac{e{\alpha }^{2}}{{ct}}{R}_{\mathrm{eq}}^{2}{B}_{\mathrm{eq}}\\ & \leqslant & (5.23\times {10}^{18}\,\mathrm{eV}){\left(\frac{{F}_{{\nu }_{p}}}{\mathrm{mJy}}\right)}^{0.84}{\left(\frac{{D}_{L}}{100\mathrm{Mpc}}\right)}^{1.67}\\ & & \times \,{\left(\frac{{\epsilon }_{e}}{{\epsilon }_{B}}\right)}^{-0.33}{f}^{-0.33}{\left(\frac{{\nu }_{p}}{\mathrm{GHz}}\right)}^{-1}{\left(\frac{t}{\mathrm{day}}\right)}^{-1}\end{array}.\end{eqnarray} \tag{ F13 }$$

The maximum energy to which protons can be accelerated during the radiatively inefficient accretion flow is limited by comparing the timescale for plasma turbulence acceleration with the diffusion timescale (Hayasaki & Yamazaki 2019). Using the photospheric radius R = R_ph, and assuming a plasma beta of 3, and a proton number density equal to the particle number density calculated above, the maximum energy of protons is (Hayasaki & Yamazaki 2019)

$$\begin{eqnarray}\begin{array}{rcl}{E}_{p} & \leqslant & \displaystyle \frac{2{\pi }^{1/2}e}{3}{\left(\displaystyle \frac{R}{{R}_{S}}\right)}^{-1}{R}_{S}{({m}_{p}{n}_{p}{c}^{2})}^{1/2}\\ {E}_{p} & \leqslant & (2.33\times {10}^{16}\,\mathrm{eV}){M}_{6}^{25/12}{m}_{\star }^{-(5/4)-(9\xi /8)}{\beta }^{2}{\left(\displaystyle \frac{t}{\mathrm{day}}\right)}^{1/6}.\end{array}\end{eqnarray} \tag{ F14 }$$