On the Diversity of Fallback Rates from Tidal Disruption Events with Accurate Stellar Structure

The tidal destruction of a star by a supermassive black hole (SMBH) generates a stream of stellar debris that, over timescales of months to years, feeds the SMBH and generates a luminous, observable signature (Hills 1975; Lacy et al. 1982; Rees 1988). These tidal disruption events (TDEs) therefore offer one of the few means to directly probe the inner regions of otherwise-quiescent galaxies, and dozens have now been observed (e.g., Gezari et al. 2012; Arcavi et al. 2014; Chornock et al. 2014; Blagorodnova et al. 2017; Hung et al. 2017; van Velzen et al. 2019; see Komossa 2015 for a review). However, our ability to confidently use the observed flares from TDEs to study, for example, black hole (BH) demographics depends critically on our physical understanding of the disruption process, and in particular the way in which the properties of the progenitor star translate to a corresponding feeding rate of the SMBH.

Along these lines, Lodato et al. (2009) simulated the full disruption of polytropes, which offer simple, yet physical descriptions of the density profiles of stellar interiors (e.g., Hansen et al. 2004), and also provided a nearly analytical means of calculating the fallback rate from a star with a given density profile (using the impulse, or "frozen-in," approximation; see also Stone et al. 2013 and Section 3 below). They found that, while the fallback rate always approached the expected, t^−5/3 scaling at late times, the peak value of the fallback and the time to peak depended on the stellar structure. Guillochon & Ramirez-Ruiz (2013) expanded this work by also studying the effect of varying the point of closest approach of the stellar center of mass to the BH, and found that denser stars more frequently leave bound "cores" that either resist the tidal shear altogether throughout pericenter passage, or reform following the full disruption of the star. The existence of these cores then modifies not only the early-time fallback, but also causes the late-time fallback to deviate from t^−5/3, and these features are a direct result of the stellar structure (in combination with the variation in the stellar pericenter).

Since then, a number of other authors have analyzed the tidal disruption of polytropes, with the aim of assessing one or another aspect of the tidal disruption process (e.g. Hayasaki et al. 2013, 2016; Coughlin & Nixon 2015; Shiokawa et al. 2015; Bonnerot et al. 2016; Coughlin et al. 2016, 2017; Bonnerot et al. 2017; Guillochon & McCourt 2017; Mainetti et al. 2017; Wu et al. 2018; Golightly et al. 2019). However, only relatively few authors have analyzed the disruption of a star with a density profile other than a polytrope. We are aware of the following studies: MacLeod et al. (2012), who investigated the disruption of giant stars by particularly massive SMBHs; Law-Smith et al. (2017), who simulated the disruption of white dwarfs with extended, hydrogen envelopes; Gallegos-Garcia et al. (2018), who analytically calculated the fallback of metal-rich material from the cores of evolved stars; and Goicovic et al. (2019), who analyzed the extent to which a more realistic stellar structure affects the stellar "disruptability."

Among the questions that remain concerning stellar structure is the degree to which more realistic stellar density profiles (i.e., those calculated with a stellar evolution code that accounts for more realistic opacities, metallicity gradients, and energy transport) affects the fallback rate compared to polytropes. More realistic stellar profiles could impose additional variability on the fallback rate, or conceivably alter characteristic timescales associated with the fallback (e.g., the time to the peak of the fallback curve). Such timescales are used to place constraints on BH properties (e.g., Guillochon & Ramirez-Ruiz 2013; Mockler et al. 2019), and hence it is necessary to understand the effects that stellar structure can have in modifying them.

In this Letter we analyze the disruption of stars with stellar structure computed with the code mesa (Paxton et al. 2011, 2013, 2015, 2018), primarily to understand the influence that such structure can have on the fallback of debris to the SMBH. In Section 2, we first describe and present the stellar profiles calculated with mesa, and in Section 3 we use the impulse approximation—as described in Lodato et al. (2009)—to calculate the fallback rate onto the BH from those stars; we show that, compared to polytropes that are matched to the same stellar mass and radius of a given mesa model, there are certain combinations of initial stellar mass and age that yield notably different fallback curves. In Section 4 we present the results of numerical simulations of disruptions of the mesa models, and compare those results to disruptions of polytropes (again, matched to the same stellar mass and radius). We discuss the implications of our simulations for estimating the BH mass from observed light curves in Section 5, and show that polytropic approximations can lead to mass-estimate discrepancies at the order of magnitude level. We summarize and conclude in Section 6.

Using the stellar evolution code mesa (Paxton et al. 2011, 2013, 2015, 2018), we evolved a 0.3 M_⊙, 1 M_⊙, and 3 M_⊙ zero-age main-sequence (ZAMS) star to the end of the main sequence (MS), being the time at which the hydrogen mass fraction in the core dropped below 0.1%. For each star, we used all of the default values for the standard inputs (e.g., each pre-MS model adopted Solar metallicity, the stars were all non-rotating, there was no mass loss in the form of winds) within mesa, version 10398.

We took snapshots of the density of each star at the ZAMS, the terminal-age main sequence (TAMS), and at one time in between ZAMS and TAMS when the hydrogen mass fraction in the core just fell below 0.2; we denote the latter by a "middle-age main sequence" star (MAMS).⁴ Figure 1 shows the density profile of each of these stars, and demonstrates, perhaps not surprisingly, that stellar evolution produces vastly different density structures over the lifetime of a given star.

Zoom In Zoom Out Reset image size

In the following two sections, we describe two different approaches to modeling the tidal disruption of these stars by an SMBH.

A useful methodology for analyzing the fallback of debris from a tidal disruption event is the impulse approximation, which posits that the tidal field of the BH acts impulsively as the center of mass of the star reaches the tidal radius. Therefore, prior to reaching the tidal radius the star retains perfect hydrostatic balance, and thereafter the star is "destroyed," meaning that each gas parcel follows its own ballistic orbit in the gravitational field of the BH. Within this approximation and to lowest order in the tidal potential, the binding energy of a given fluid parcel is only a function of the projected distance of that fluid parcel from the BH onto the line connecting the stellar center of mass and the BH. This energy dependence implies that perpendicular "slices" of the star return simultaneously to the BH, which allows one to construct a fallback rate $\dot{M}$—being the rate at which bound material returns to the BH—that accounts for the stellar density profile ρ; the result is (Lodato et al. 2009; Gallegos-Garcia et al. 2018; Golightly et al. 2019)

$$\begin{eqnarray}&&\dot{M}=\displaystyle \frac{{M}_{\star }}{{T}_{\mathrm{mb}}}{\tau }^{-5/3}{\int }_{{\tau }^{-2/3}}^{1}\displaystyle \frac{\rho (x){xdx}}{{\rho }_{\star }},\end{eqnarray} \tag{ 1 }$$

where M_⋆ is the stellar mass, τ = t/T_mb with

$$\begin{eqnarray}&&{T}_{\mathrm{mb}}={\left(\displaystyle \frac{{R}_{\star }}{2}\right)}^{3/2}\displaystyle \frac{2\pi M}{{M}_{\star }\sqrt{{GM}}}\end{eqnarray} \tag{ 2 }$$

the return time of the most bound debris, M the BH mass and R_⋆ the stellar radius, η = R/R_⋆ with R the spherical distance from the center of the star, and ${\rho }_{\star }=3{M}_{\star }/(4\pi {R}_{\star }^{3})$ is the average stellar density. Note that the integral in Equation (1) is negative for all times t < T_mb, and hence the physical fallback rate is zero until the most bound debris element returns to the BH.

From Equation (1), within the impulse approximation the only timescale relevant to the problem is the return time of the most bound debris; hence the time at which the fallback reaches any characteristic value (e.g., the time to peak, the time between half-peaks, the time to reach $\dot{M}\propto {t}^{-5/3+\alpha }$ with α > 0) is also a constant multiple—which depends on the stellar structure—of this timescale. If we denote the dimensionless time to any characteristic fallback value by τ_c, then by definition the corresponding physical time t_c is

$$\begin{eqnarray}&&{t}_{{\rm{c}}}={\left(\displaystyle \frac{{R}_{\star }}{2}\right)}^{3/2}\displaystyle \frac{2\pi M}{{M}_{\star }\sqrt{{GM}}}{\tau }_{{\rm{c}}}(\rho ).\end{eqnarray} \tag{ 3 }$$

This simple expression demonstrates the relative importance that stellar structure can have when inferring BH properties from observations⁵ : assume that for a given tidal disruption event with known t_c, we know the mass and radius of the disrupted star. Then, for the same event if we ascribe to the disrupted star two different density profiles, ρ₁ and ρ₂, we will infer two different BH masses for the same event, M₁ and M₂, their ratio being

$$\begin{eqnarray}&&\displaystyle \frac{{M}_{2}}{{M}_{1}}={\left(\displaystyle \frac{{\tau }_{{\rm{c}}}({\rho }_{1})}{{\tau }_{{\rm{c}}}({\rho }_{2})}\right)}^{2}.\end{eqnarray} \tag{ 4 }$$

Thus, relatively small differences in the stellar structure are capable of producing more discrepant BH masses owing to the squared dependence in this expression.

Figure 2 illustrates the fallback rate onto the BH (in units of Solar masses per year as a function of time in years) from the frozen-in approximation, where the disrupted star has a mass and radius equal to those of the 1 M_⊙, TAMS star (see Table 1) and the BH has a mass of 10⁶ M_⊙. Each curve adopts a different stellar density profile, being a γ = 5/3 polytrope (blue, dotted–dashed), a γ = 1.35 polytrope (green, dashed), and the profile resulting from mesa (red, solid). It is clear from this figure that, despite the fact that the bulk stellar properties are the same, the density profile has a marked effect on the shape of the curve. To highlight the differences induced by stellar structure, the points on each curve give the time to different, characteristic values of the fallback curve, being the time to peak (t_max), the time to reach half the peak (which occurs on the rise and the decay of the curve, denoted, respectively, by t_half,1 and t_half,2), and the time to reach $\dot{M}\propto {t}^{-4/3}$, t_4/3. While certain characteristic times are visibly comparable for each model (e.g., t_max), other timescales are more noticeably discrepant (e.g., t_4/3).

Zoom In Zoom Out Reset image size

Because of these discrepancies, for a TDE with a given physical timescale, the BH mass required to yield a fallback curve with that timescale will differ according to the stellar structure model that one employs, and the magnitude of the difference in mass will be larger or smaller depending on the characteristic timescale itself. For example, using Equation (4) and letting ρ₁ describe the density profile of the 1 M_⊙, TAMS mesa progenitor, we find ${M}_{2}/{M}_{1}\simeq 2.3$ if one models the density profile by a γ = 1.35 polytrope and uses the difference between the time to peak and the first time to half peak (i.e., the timescale t_c = t_max–t_half,1⁶ ); physically, one can interpret this result by saying that a γ = 1.35 polytrope requires a BH mass ∼2.3 times larger to reproduce the time to peak from the 1 M_⊙, TAMS mesa progenitor. On the other hand, if we use the same timescale but model the star as a γ = 5/3 polytrope, then we find ${M}_{2}/{M}_{1}\simeq 14$—because a γ = 5/3 polytrope peaks considerably earlier than the real star, a much larger BH mass is required in order to extend the time to maximum fallback.

Figure 3 directly demonstrates how changing the mass of the BH yields commensurate characteristic timescales: the solid red curve shows the same fallback rate as in Figure 2 for the 1 M_⊙, TAMS stellar profile computed with mesa and a 10⁶ M_⊙ SMBH (note that this figure is on a linear–linear scale). The dashed, green line, and the dotted–dashed blue line show the fallback rates for a γ = 1.35 and γ = 5/3 polytrope, respectively, again with the same stellar properties as those of the mesa progenitor. However, here we varied the mass of the BH according to the value required to yield the same time to peak, being M = 2.3 × 10⁶ for the γ = 1.35 polytrope and $M=1.4\times {10}^{7}{M}_{\odot }$ for the γ = 5/3 polytrope. The time between the first half-peak (marked with a †) and the peak fallback rate (marked with a ⋆) is the same in each case. This shows that, by varying the BH mass appropriately, one can reconstruct the same physical timescale for different stellar properties.

Zoom In Zoom Out Reset image size

One can repeat this procedure for the nine different models presented in Table 1 and thereby assess the degree to which a polytropic density profile reproduces the fallback curve—and correspondingly the inferred BH mass—of the one obtained with the mesa model. However, the impulse approximation ignores crucial physics of the disruption process that also alter the fallback rate (Guillochon & Ramirez-Ruiz 2013; Coughlin & Nixon 2015; Steinberg et al. 2019). In the next section, we employ hydrodynamical simulations to obtain more realistic fallback rates.

4.1. Setup

4.2. Results

Figure 4 illustrates the fallback rate onto the BH in Solar masses per year as a function of time in years from six different stellar models, with the specific stellar model shown in the legend. In each panel the solid curve shows the fallback from the star with the mesa density profile, the dashed curve is the γ = 5/3 model matched to the mesa stellar mass and radius, and the dotted–dashed curve is the γ = 1.35 polytropic model (again, with the same mass and radius as the mesa star). It is evident from the top-left panel and the middle-left panel that the 0.3 M_⊙, ZAMS and the 1.0 M_⊙, ZAMS mesa fallback curves are extremely well-reproduced by γ = 5/3 and γ = 1.35 polytropes, respectively. This finding indicates, correspondingly, that such stars can be very well-modeled by single polytropes that are gas-pressure and (nearly) radiation-pressure dominated (see Figure 5). This result for the 1 M_⊙, ZAMS star was also recovered by Goicovic et al. (2019), who found that their fallback rates (from a 1 M_⊙, ZAMS star generated with mesa) were very similar to those of Guillochon & Ramirez-Ruiz (2013), who used a polytropic approximation.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

For every other stellar model, however, there are notable differences between the fallback curves from the mesa and the polytropic models. Specifically, we see that employing the mesa density profile over either polytropic model systematically shifts the return time of the most bound debris to earlier times, the time to peak to earlier times, and the magnitude of the peak rate itself is also increased. Furthermore, because the total mass is the same, the larger peak in the accretion rate and the earlier time-to-peak imply that the mesa curves must fall below the polytropic ones at some point, and this is indeed recovered in each case. It is also apparent that the mesa models conform to a power-law decline at an earlier time than do the polytropic models.

Every polytropic star is completely disrupted by the SMBH for these β = 3 encounters, which agrees with the results of Guillochon & Ramirez-Ruiz (2013) and Mainetti et al. (2017), who demonstrated that the critical β for the full disruption of a γ = 5/3 polytrope is $\beta \simeq 0.9$, while that for a γ = 4/3 polytrope is $\beta \simeq 2$. Each ZAMS mesa progenitor is also fully disrupted. However, more highly evolved stars start to yield partial TDEs, and the mesa 1 M_⊙ MAMS and 3.0 M_⊙ MAMS leave stellar cores at the location of the marginally bound orbit. We see that, for the case of the 3.0 M_⊙, MAMS progenitor, the presence of the core has the affect of modifying the late-time fallback rate, which declines approximately as ∝ t^−9/4 instead of t^−5/3. This result is in agreement with the analytic predictions of Coughlin & Nixon (2019).

The TAMS mesa progenitors all possess extremely dense cores that survive the tidal encounter, each of which modifies the late-time fallback rate onto the BH, as shown in Figure 6 (the time taken for the 3.0 M_⊙ progenitor to go from MAMS to TAMS is very short, and hence the TAMS fallback curve appears nearly identical to the bottom-right panel of Figure 4; we therefore opted not to show this fallback rate). We also emphasize that the lifetime of the 0.3 M_⊙ star is well in excess of the age of the universe, and hence this star cannot be disrupted by an SMBH (at least not any time soon). However, the density profile of this star could conceivably be achieved by a more massive progenitor, which would evolve to the TAMS within the age of the universe, with different initial properties (e.g., metallicity, rotation).

Zoom In Zoom Out Reset image size

As we noted above, the 1 M_⊙ MAMS mesa progenitor also possesses a bound stellar core at the center of mass of the tidally disrupted debris stream. In this case, however, the star is initially completely disrupted, and the core recollapses out of the stream at a time significantly after the stellar center of mass passes through pericenter. For this reason, the mass of the surviving core is only a small fraction of the initial progenitor ($\simeq 15 \%$), and consequently the fallback curve shows little evidence of the gravitational influence of the core over ∼1 yr and appears to approach a decline ∝ t^−5/3.

Figure 7 shows the fallback from the 1 M_⊙ MAMS mesa progenitor out to 10 yr post-disruption and demonstrates, however, that the presence of the core does start to affect the fallback at later times. In particular, we see that while the first year shows little evidence of the existence of a bound core, there is a clear break in the fallback curve at a time of 1–2 yr where the power law of the fallback rate transitions from ∝ t^−5/3 to one that is better matched by ∝ t^−9/4. Interestingly, this time at which a break in the power law is exhibited is very close to the time predicted by the analytical model in Coughlin & Nixon (2019; see their Figure 2), and coincidentally also occurs around the same time at which the fallback rate drops below the Eddington limit of the SMBH (dashed black line, assuming a radiative efficiency of 10% and an electron-scattering opacity of 0.34 cm² g⁻¹).

Zoom In Zoom Out Reset image size

It is apparent from Figures 4 and 6 that, for the majority of progenitors, non-polytropic stellar structure generates substantial differences in the fallback rate onto the BH. Notably, the mesa profiles yield earlier times to peak and larger peak fallback rates, and they more rapidly approach a power-law falloff as compared to the polytropes. To use these differences to estimate the corresponding differences in the inferred BH mass that would arise by assuming a given density profile, we must have a mapping between a characteristic timescale (e.g., the time to peak), the BH mass, and the properties of the star. As shown in Section 3, when the fallback rate is computed with the impulse approximation, this mapping arises through the timescale

$$\begin{eqnarray}&&{T}_{\mathrm{mb}}\simeq \displaystyle \frac{2\pi {R}_{\star }^{3/2}M}{{M}_{\star }\sqrt{{GM}}}\propto {M}^{1/2},\end{eqnarray} \tag{ 5 }$$

which is the return time of the most bound debris. Additional dependence on the stellar structure modifies the fallback rate through a dimensionless function of time normalized by T_mb, and that dimensionless function can be calculated from the (assumed-unaltered) density profile when the star is at the tidal radius. Thus, the dependence on the BH mass arises only as ∝ M^1/2, and in this approximation, fallback curves are simply scaled in time and magnitude by $\sqrt{M}$.

By comparing Figure 2 to Figures 4 and 6, wee see that the frozen-in approximation does not accurately reproduce many of the features of the numerically obtained fallback rates. In addition to the fact that the time to peak is shorter and the peak itself is higher in the numerical simulations (by a factor ≳10), the ordering of the curves is actually inverted between the two approaches: while Figure 2 shows that the γ = 5/3 polytrope peaks earlier than γ = 1.35 polytrope, which itself peaks earlier than the mesa model, Figures 4 and 6 demonstrate that the γ = 5/3 polytrope always reaches a peak after the γ = 1.35 polytrope.⁷ Moreover, for every case except the 0.3 M_⊙ ZAMS and 0.3 M_⊙ MAMS progenitors, the mesa model peaks earlier than the γ = 1.35 polytrope. The numerically obtained return time of the most bound debris also differs for each density profile, whereas, under the frozen-in approximation, this timescale—for the same BH mass—is only affected by the stellar mass and radius (which, for a given star, are identical by construction).

These discrepancies indicate that the impulse approximation does not include enough physics to accurately capture the bulk features of the fallback rate. As discussed at length in Coughlin et al. (2016), it is likely that the most crucial physical ingredient lacking from the impulse approximation is the self-gravity of the debris stream, as the stellar center of mass rapidly recedes outside of the tidal sphere of the BH. At this point, the stellar density is comparable to the "BH density," being ${\rho }_{\bullet }\sim M/{r}_{{\rm{t}}}^{3}$, and the self-gravity of the stream is capable of competing against the shear of the BH. The self-gravity of the stream induces density waves that traverse the stream radially, and these waves serve to generate more pronounced "shoulders" near the extremities of the stream and flatten the dm/dr ∝ ρ curve from the polytropic one that follows from the frozen-in approximation. It is, in fact, because of this nearly flat dm/dr generated by self-gravity that the fallback curves more rapidly approach the t^−5/3 decline (or the t^−9/4 decline for the partial TDEs). The higher central density of the γ = 1.35 polytrope also generates more vigorous density waves, which correspondingly produce a flatter density distribution and give rise to an earlier time-to-peak as compared to a γ = 5/3 polytrope.

Nonetheless, it is likely that after some amount of time following the disruption, the mass distribution is approximately frozen-in, meaning that self-gravity has smoothed out any density perturbations and the stream is long enough that the time-dependent potential due to self-gravity is small.⁸ In this case, the energies of gas parcels comprising the stream are still Keplerian in the potential of the BH, and the energies themselves are simply established at some later time; indeed, these arguments were used and verified by a direct evaluation of the energy distribution at different times by Lodato et al. (2009) and Guillochon & Ramirez-Ruiz (2013) to calculate fallback rates to the BH after only a small fraction of the return time of the most bound debris had been directly simulated.⁹ Moreover, if for a given β and stellar progenitor the density profile at the time the energy is frozen-in is independent of the BH mass, which the simulations of Wu et al. (2018) verify¹⁰ (see their Figure 1), then it follows that any physical timescale in a TDE can be written

$$\begin{eqnarray}&&{t}_{{\rm{c}}}=\sqrt{M}{f}_{\star ,\beta },\end{eqnarray} \tag{ 6 }$$

where f_⋆,β is a function that depends only on the stellar properties and β. We therefore see that, for the same star and the same orbital parameters, we recover the same result as we did in Section 3: for two TDEs with identical orbital and stellar properties and physical timescales t₁ and t₂, we can satisfy t₁ = t₂ by changing the BH mass M₁ to M₂, with M₂ given by

$$\begin{eqnarray}&&{M}_{2}={M}_{1}{\left(\displaystyle \frac{{t}_{2}}{{t}_{1}}\right)}^{2}.\end{eqnarray} \tag{ 7 }$$

As an example, the 3 M_⊙, MAMS progenitor (bottom-right panel of Figure 4) modeled as a γ = 1.35 polytrope has a time to go from the first-half-max, t_half,1, to max, t_max, of ${t}_{1}\,={t}_{\max }\mbox{--}{t}_{\mathrm{half},1}\simeq 0.085\,\mathrm{yr}$. On the other hand, the mesa model of the same star has ${t}_{2}={t}_{\max }\mbox{--}{t}_{\mathrm{half},1}=0.0275\,\mathrm{yr}$. Thus, if we modeled the disruption of the mesa star as a polytrope, then we would require a BH mass of ${M}_{2}={({t}_{2}/{t}_{1})}^{2}{M}_{1}\,\simeq 0.097{M}_{1}\,\simeq {10}^{5}{M}_{\odot }$ to reproduce the observed timescale.

Table 2 gives the ratio M₂/M₁ required to shift the timescale of the polytropic star to the timescale reproduced by the mesa-star disruption. The timescale itself is shown in the top row of the table, where t_max–t_half,1 is the time to go from first-half-max to the peak fallback rate, t_half,2–t_max is the time taken to fall by a factor of two below the peak fallback rate, and t_half,2–t_half,1 is the FWHM of the fallback curve. The stellar progenitor is given in the left-most column of the table, and the value in the left (right) of each cell is the ratio M₂/M₁ required to yield the mesa-generated timescale by modeling the star as a γ = 1.35 (γ = 5/3) polytrope. For example, if we were to model the 0.3 M_⊙ star as a γ = 1.35 polytrope, then we would require a BH mass of M₂ = 5 × M₁ to reproduce the timescale t_max–t_half,1 found from the disruption of the mesa model, and the number M₂/M₁ = 5 is shown in the top-left cell of the table.

Table 2.  The Ratio M₂/M₁ that is Required to Produce the Same Physical Timescale if the mesa Fallback Curve is Modeled by a Polytropic One

Timescale	${t}_{\max }\mbox{--}{t}_{\mathrm{half},1}$		${t}_{\mathrm{half},2}\mbox{--}{t}_{\max }$		${t}_{\mathrm{half},2}\mbox{--}{t}_{\mathrm{half},1}$
Star	γ = 1.35	γ = 5/3	γ = 1.35	γ = 5/3	γ = 1.35	γ = 5/3
0.3 M⊙ ZAMS	M2/M1 = 5.0	M2/M1 = 0.81	M2/M1 = 12	M2/M1 = 1.19	M2/M1 = 9.1	M2/M1 = 1.07
0.3 M⊙ MAMS	M2/M1 = 7.0	M2/M1 = 1.2	M2/M1 = 5.3	M2/M1 = 0.50	M2/M1 = 5.9	M2/M1 = 0.67
0.3 M⊙ TAMS	M2/M1 = 0.24	M2/M1 = 0.061	M2/M1 = 0.22	M2/M1 = 0.022	M2/M1 = 0.23	M2/M1 = 0.030
1.0 M⊙ ZAMS	M2/M1 = 0.83	M2/M1 = 0.15	M2/M1 = 0.92	M2/M1 = 0.074	M2/M1 = 0.89	M2/M1 = 0.093
1.0 M⊙ MAMS	M2/M1 = 0.33	M2/M1 = 0.042	M2/M1 = 0.25	M2/M1 = 0.025	M2/M1 = 0.27	M2/M1 = 0.031
1.0 M⊙ TAMS	M2/M1 = 0.14	M2/M1 = 0.032	M2/M1 = 0.13	M2/M1 = 0.010	M2/M1 = 0.14	M2/M1 = 0.015
3.0 M⊙ ZAMS	M2/M1 = 0.19	M2/M1 = 0.033	M2/M1 = 0.22	M2/M1 = 0.019	M2/M1 = 0.21	M2/M1 = 0.023
3.0 M⊙ MAMS	M2/M1 = 0.097	M2/M1 = 0.018	M2/M1 = 0.059	M2/M1 = 0.0057	M2/M1 = 0.072	M2/M1 = 0.0086
3.0 M⊙ TAMS	M2/M1 = 0.067	M2/M1 = 0.0076	M2/M1 = 0.045	M2/M1 = 0.0049	M2/M1 = 0.053	M2/M1 = 0.0058

Note. Here the physical timescale used to infer the required mass ratio is given in the top row, with T_max–t_half,1 the time from the first-half-max to the maximum fallback rate, t_half,2–t_max the time taken to go from the peak fallback rate to half that value, and t_half,2–T_half,1 the FWHM; each one of these timescales differs for a given star, but by scaling the BH mass by the value shown in each cell, they can be brought into agreement with one another. The stellar model is given in the left-most column, and the ratio M₂/M₁ obtained by using a γ = 1.35 (γ = 5/3) is given in the left (right) of each cell.

Download table as:  ASCIITypeset image

We see from this table that, for stellar density profiles that are well-reproduced by polytropes (the 0.3 M_⊙ ZAMS and the 1 M_⊙ ZAMS stars, as shown in Figure 5), the mass ratio required to reproduce the mesa-generated timescale with a polytropic one is very close to unity. In these instances, the inferred BH mass estimate is not far off the true, underlying value. However, for more massive progenitors and more highly evolved stars, the extremely dense core of the mesa progenitor shifts each timescale earlier, which consequently requires a significantly smaller BH mass to yield the same characteristic timescale with a polytropic model.

As a direct demonstration of this effect, Figure 8 illustrates the fallback rate from the disruption of the 3 M_⊙, ZAMS mesa model by an SMBH of mass M₁ = 10⁶M_⊙, which is the same curve shown in the bottom-left panel of Figure 4. The dotted–dashed curve is the fallback curve from the disruption of a γ = 1.35 polytrope, the stellar mass and radius identical to those of the mesa star. In this case, however, we set the mass of the disrupting SMBH to M₂ = 2 × 10⁵, which is, as seen in Table 2, the value predicted to equate the time to go from the first-half-max to the peak between the two models, and we aligned the first time to half-max of the polytrope fallback curve to that of the mesa star (i.e., we initially "see" both disruptions at the same time, that time being the first time to half-max). We see that this polytropic model provides an extremely good fit to the "data" obtained from the fallback curve of the mesa model, but at the expense of incorrectly inferring the SMBH mass by nearly an order of magnitude.

Zoom In Zoom Out Reset image size

Of course, our approach here to "modeling" the light curve of a TDE is overly simplistic, as one does not necessarily have any prior information about the nature of the progenitor or the BH, and one must use a combination of timescales to recover the best-fit model parameters (as is done in, for example, Guillochon et al. 2018). However, these results do demonstrate that care needs to be taken to ensure that the template fallback curves used to interpret observed data sets contain a sufficiently broad range of stellar density profiles (e.g., accounting for stellar age, and non-polytropic density profiles) to ensure that the inferred parameters are accurate and that the error bars that result from the data fitting are appropriate.

In this Letter we presented simulations of the tidal disruption of stars encountering an SMBH. We modeled the tidally disrupted stars with the same bulk properties (mass and radius) using three different prescriptions: (1) a γ = 5/3, polytropic density profile, (2) a γ = 1.35(≈4/3), polytropic density profile, and (3) a density profile calculated from the mesa stellar evolution code. Our simulated disruptions included stars with masses of 0.3 M_⊙, 1.0 M_⊙, and 3.0 M_⊙, each of which was evolved to the ZAMS, the TAMS (where the hydrogen mass fraction in the core fell below 0.1%), and the "middle-age main sequence," which we defined to be the time at which the hydrogen mass fraction in the core fell below 20%. We therefore simulated a total of 27 disruptions (nine different stars, each star modeled with a mesa density profile and two different polytropic profiles).

In each simulation we maintained the same physics: we employed a polytropic equation of state where the Lagrangian entropy was fixed and set to ensure the isolated star was in hydrostatic balance (see Figure 9); we fixed the adiabatic index in each simulation to 5/3, which ensures that the dynamics of the stream's self-gravity (see Coughlin & Nixon 2015) differs only by the mass-entropy distribution along the stream; and we fixed the tidal effects from the BH on each star by employing a $\beta \equiv {r}_{{\rm{t}}}/{r}_{{\rm{p}}}=3$ for each simulation, where r_p (r_t) is the pericenter (tidal) radius of the star (see also Table 1). Our simulations therefore isolate the impact of the stellar density profile calculated from mesa when compared to those calculated by γ = 5/3 and γ = 1.35 polytropes.

In general we find that there are significant differences in the simulated fallback rates for stars with different masses and different ages, and further that in most cases these fallback rates deviate significantly from predictions made using polytropes. The exceptions, which come as no surprise as their structures are accurately modeled by polytropes (see Figure 5), are the 0.3 M_⊙ ZAMS star—which is well-modeled by a γ = 5/3 polytrope—and the 1.0 M_⊙ ZAMS star—which is well-modeled by a γ = 1.35 polytrope. At both MAMS and TAMS we find that neither polytrope provides an acceptable description for the fallback curve. Similarly the fallback rates from the 3.0 M_⊙ stars are all significantly different to the fallback rates from either polytrope.

There are also differences found in the overall dynamics of the disruption event. In several cases, most notably for the 0.3 M_⊙ MAMS star, the debris stream is significantly more self-gravitating for the mesa star than for the polytropes. This results in the fallback curve exhibiting more variability on the power-law decay (see the top-right panel of Figure 4). We also find that several of the mesa stars are not fully disrupted, even with an impact parameter of β = 3, and this arises from the more centrally concentrated nature of the mesa stars.¹¹ The difficulty of fully disrupting real stars, and particularly more highly evolved stars, implies that a greater fraction of the events we observe will be partial, rather than full, disruptions (though we note that stars spend more of their lives near the ZAMS, where the mesa profiles still yield full disruptions for β = 3). It has recently been shown (Coughlin & Nixon 2019) that TDEs that leave a bound core have a fallback rate whose power-law index asymptotes to ≈−9/4 rather than the usual −5/3. In each case that leaves a bound core, our simulations recover this result. In future work we will explore whether simulations are also capable of recovering, e.g., the time at which the power-law slope changes to this value as a function of the mass of the core that survives the encounter. We find (see also Guillochon & Ramirez-Ruiz 2013) that in some cases, here for the 1.0 M_⊙ MAMS mesa star, that the core can be initially fully disrupted but reform after leaving the tidal radius. We attribute this to the velocity field imparted in the stellar debris by the BH tides, which can at later times cause the stream to converge along its width and augment the density within the stream (Coughlin et al. 2016; Steinberg et al. 2019).

In Section 5 we described the impact of using realistic stellar models in simulations of TDEs on the inference of the BH mass from observed TDE light curves. We showed that for many types of stars, and particularly those that are more massive at the ZAMS or more highly evolved, the shape of the fallback curves from simulations that employ polytropic stellar models can lead to large errors in the estimated BH mass. Therefore, employing accurate models for the imprint of stellar structure on the fallback rate to produce templates for TDE light curves appears essential for accurately inferring the BH mass. It could also be that degeneracies between the parameters in TDEs (e.g., stellar mass, stellar age, stellar spin, stellar metallicity, BH mass, BH binarity, BH spin, impact parameter, inclination and phase angles of the stellar orbit, and the accretion dynamics on small scales) mean that an unambiguous estimation of system parameters (or at least an understanding of the true level of error within those estimates) requires accurate and detailed modeling of the stars that are being disrupted.

We thank the anonymous referee for a useful report. C.J.N. is supported by the Science and Technology Facilities Council (grant No. ST/M005917/1). E.R.C. acknowledges support from the Lyman Spitzer Jr. Fellowship and NASA through the Einstein Fellowship Program, grant PF6-170170, and the Hubble Fellowship, grant No. HST-HF2-51433.001-A awarded by the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Incorporated, under NASA contract NAS5-26555. This work was performed using the DiRAC Data Intensive service at Leicester, operated by the University of Leicester IT Services, which forms part of the STFC DiRAC HPC Facility (www.dirac.ac.uk). The equipment was funded by BEIS capital funding via STFC capital grants ST/K000373/1 and ST/R002363/1 and STFC DiRAC Operations grant ST/R001014/1. DiRAC is part of the National e-Infrastructure.

Software: splash (Price 2007), phantom (Price et al. 2018), mesa (Paxton et al. 2011, 2013, 2015, 2018).

Figure 9 illustrates a comparison between the density profiles obtained with mesa (red curves) and those obtained with phantom after the initial particle distribution has been relaxed (black curves). For the 1 M_⊙ TAMS star, the phantom profile overshoots the mesa one by about 10%, and the density of the very outermost radii of the 3 M_⊙ star (at all ages) is slightly larger. However, in most cases the two curves are nearly indistinguishable over the entire range in radius of the star.

Zoom In Zoom Out Reset image size