# The fine line between total and partial tidal disruption events

Deborah Mainetti, Alessandro Lupi, Sergio Campana, Monica Colpi, Eric R.

Coughlin, James Guillochon, Enrico Ramirez-Ruiz

► To cite this version:

Deborah Mainetti, Alessandro Lupi, Sergio Campana, Monica Colpi, Eric R. Coughlin, et al. The fine line between total and partial tidal disruption events. Astronomy & Astrophysics - A&A, 2017, 600, pp.A124. 10.1051/0004-6361/201630092. hal-01554841

HAL Id: hal-01554841

https://hal.science/hal-01554841v1

Submitted on 26 Jun 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# The fine line between total and partial tidal disruption events

Deborah Mainetti<fnref n="1" />,<fnref n="2" />,<fnref n="3" />, Alessandro Lupi<fnref n="3" />,<fnref n="4" />, Sergio Campana<fnref n="2" />, Monica Colpi<fnref n="1" />,<fnref n="3" />, Eric R. Coughlin<fnref n="5" />,*, James Guillochon<fnref n="6" />, and Enrico Ramirez-Ruiz<fnref n="7" />

Received 18 November 2016 / Accepted 24 February 2017

## ABSTRACT

Flares from tidal disruption events are unique tracers of quiescent black holes at the centre of galaxies. The appearance of these flares is very sensitive to whether the star is totally or partially disrupted, and in this paper we seek to identify the critical distance of the star from the black hole ($r_d$) that enables us to distinguish between these two outcomes. We perform here mesh-free finite mass, traditional, and modern smoothed particle hydrodynamical simulations of star-black hole close encounters, with the aim of checking if the value of $r_d$ depends on the simulation technique. We find that the critical distance (or the so-called critical disruption parameter $\beta_d$) depends only weakly on the adopted simulation method, being $\beta_d = 0.92 \pm 0.02$ for a $\gamma = 5/3$ polytrope and $\beta_d = 2.01 \pm 0.01$ for a $\gamma = 4/3$ polytrope.

**Key words.** hydrodynamics – methods: numerical – galaxies: nuclei – black hole physics – accretion, accretion disks

## 1. Introduction

There is compelling evidence of the ubiquitous presence of massive black holes (BHs) at the centres of nearby galaxies (Kormendy & Richstone 1995; Kormendy & Ho 2013). These BHs are mainly in low-luminous states (Ho 2008) or in quiescence, but sometimes they can enter highly luminous phases (AGN) that are due to sudden influxes of the surrounding gas. These influxes can be provided by the tidal disruption (TD) of stars (Rees 1988). TDs occur when stellar dynamical encounters scatter a star (of mass $M_*$ and radius $R_*$) onto a low angular momentum orbit about the BH (of mass $M_{BH}$), subjecting it to the extreme BH tidal field (Alexander 2012). Specifically, if the star comes close to the so-called BH tidal radius

$$r_t \sim R_* \left( \frac{M_{BH}}{M_*} \right)^{1/3} \sim 10^2 R_* \left( \frac{R_*}{1 R_*} \right) \left( \frac{M_{BH}}{10^6 M_*} \right)^{1/3} \left( \frac{1 M_*}{M_*} \right)^{1/3} \quad (1)$$

(Hills 1975; Frank & Rees 1976), the star will be totally or partially disrupted, depositing a fraction of its mass onto the BH through an accretion disc and powering a bright flare (e.g. Rees 1988; Phinney 1989; Evans & Kochanek 1989; Lodato et al. 2009; Strubbe & Quataert 2009; Lodato & Rossi 2011; Guillochon & Ramirez-Ruiz 2013, 2015a; Coughlin & Begelman 2014). For a star to be disrupted outside the event horizon

of a BH, that is, in order to observe the corresponding TD accretion flare, $r_t$ must be greater than the BH event horizon radius

$$r_s = \frac{xGM_{BH}}{c^2} \sim 4 R_* \left( \frac{M_{BH}}{10^6 M_*} \right) \left( \frac{x}{2} \right), \quad (2)$$

where $x$ encapsulates effects related to the BH spin (Kesden 2012). Hence, the non-rotating destroyer BH mass must be $M_{BH} \lesssim 10^8 M_*$ when solar-type stars are involved. TD accretion flares thus reveal otherwise quiescent or low-luminous BHs in a mass range complementary to that probed in AGN surveys (Vestergaard & Osmer 2009).

Regardless of whether the destruction is total or partial, most of the stars in a galaxy fated to be disrupted by the central BH are scattered onto low angular momentum orbits from about the BH sphere of influence, that is, onto nearly parabolic trajectories (Magorrian & Tremaine 1999; Wang & Merritt 2004). For this reason, in this work we assume the disruptive orbits to be parabolic. This assumption, together with the kick naturally imparted by the disruption itself (Manukian et al. 2013), prevents our partially disrupted stars from encountering the BH a multitude of times. In this parabolic regime, about half of the stellar debris produced by a (total or partial) stellar disruption binds to the BH, returns to pericentre on different elliptical orbits (that is, with different orbital energies; Lacy et al. 1982), circularises forming an accretion disc, and falls back onto the

<footnote n="1">¹ Dipartimento di Fisica G. Occhialini, Università degli Studi di Milano Bicocca, Piazza della Scienza 3, 20126 Milano, Italy
e-mail: d.mainetti1@campus.unimib.it</footnote>

<footnote n="2">² INAF, Osservatorio Astronomico di Brera, via E. Bianchi 46, 23807 Merate (LC), Italy</footnote>

<footnote n="3">³ INFN, Sezione di Milano-Bicocca, Piazza della Scienza 3, 20126 Milano, Italy</footnote>

<footnote n="4">⁴ Institut d'Astrophysique de Paris, Sorbonne Universités, UPMC Univ Paris 6 et CNRS, UMR 7095, 98 bis bd Arago, 75014 Paris, France</footnote>

<footnote n="5">⁵ Astronomy Department and Theoretical Astrophysics Center, University of California, Berkeley, CA 94720, USA</footnote>

<footnote n="6">⁶ Harvard-Smithsonian Center for Astrophysics, The Institute for Theory and Computation, 60 Garden Street, Cambridge, MA-02138, USA</footnote>

<footnote n="7">⁷ Department of Astronomy & Astrophysics, University of California, Santa Cruz, CA-95064, USA</footnote>

<footnote>* Einstein Fellow.</footnote>
BH emitting a peculiar flare. The fallback rate is likely to be somewhat different from the rate of debris returning to pericentre (e.g. Cannizzo et al. 1990; Ramirez-Ruiz & Rosswog 2009; Hayasaki et al. 2013, 2016; Coughlin & Nixon 2015; Guillochon & Ramirez-Ruiz 2015b; Piran et al. 2015; Shiokawa et al. 2015; Bonnerot et al. 2016; Coughlin et al. 2016), which in turn depends on the structure of the disrupted star (e.g. Lodato et al. 2009) and the properties of the encounter (e.g. Guillochon & Ramirez-Ruiz 2013, 2015a).

In this paper we aim at computing the critical distance $r_d$ at which a star becomes totally or partially disrupted by the BH tidal field. We are interested in finding the critical disruption parameter

$$ \beta_d = \frac{r_t}{r_d} = \beta \frac{r_p}{r_d} \qquad (3) $$

for specific stellar structures that distinguishes total TDs from partial TDs, with $r_t$ given by Eq. (1), $\beta = r_t/r_p$ and $r_p$ being the pericentre of the star around the BH. A partial TD is obtained for $\beta < \beta_d$, that is, for $r_p > r_d$, a total TD for $\beta \ge \beta_d$, that is, for $r_p \le r_d$. The need to introduce the critical distance $r_d$ arises because the tidal radius $r_t$ defines where the BH tidal force overcomes the stellar self-gravity only at the stellar surface, and not everywhere within the star. This problem has been considered previously (Guillochon & Ramirez-Ruiz 2013, 2015a; hereafter GRR). GRR evaluated the critical disruption parameter $\beta_d$ for polytropes of index 5/3 and 4/3 (which represent low- and high-mass stars, respectively) using a series of adaptive mesh refinement (AMR) grid-based hydrodynamical simulations of tidal encounters of star and BH. In this paper, we instead present the results of simulations we performed for the same purpose with the codes GADGET2 (traditional smoothed particle hydrodynamical (SPH); Springel 2005)<fnref n="1" /> and GIZMO (modern SPH and mesh-free finite mass (MFM); Hopkins 2015)<fnref n="2" />. Since these techniques all have advantages but also limits, we are inclined to follow GRR in finding the critical disruption parameter $\beta_d$ for certain stellar structures<fnref n="3" /> using an MFM, a traditional SPH, and a modern SPH code instead of an AMR grid-based code, with the goal of comparing results from different techniques.

The paper is organised as follows. In Sect. 2 we compare AMR grid-based codes to GIZMO MFM, traditional SPH, and modern SPH techniques. In Sect. 3 we discuss our method and describe how we evaluate the stellar mass loss $\Delta M$ in our simulated encounters. We show the curves of mass loss we obtain for all codes as a function of $\beta$ and polytropic index, comparing them and the corresponding $\beta_d$ with GRR. Section 4 contains our summary.

## 2. Grid-based vs. SPH and GIZMO MFM codes

Fluid hydrodynamics and interactions in astrophysics are generally treated using two different classes of numerical methods: Eulerian grid-based (e.g. Laney 1998; Leveque 1998) and Lagrangian SPH (e.g. Monaghan 1992; Price 2005, 2012; Cossins 2010). Basically, grid-based methods divide a domain into stationary cells traversed over time by the investigated fluid, and account for information exchange between adjacent cells in

the aim at solving the fluid equations. In particular, AMR grid-based techniques (e.g. Berger & Oliger 1984; Berger & Colella 1989) adapt the cell number and size according to the properties of different fluid regions, thus increasing the resolution where needed (for example in high-density regions) and reducing computational efforts and memory employment where lower resolution is sufficient. In contrast, SPH methods are Lagrangian by construction and model a fluid as a set of interacting fluid elements, or particles, each with its own set of fluid properties. In practice, the density of each particle is calculated by considering the neighbours within its so-called smoothing length (e.g. Price 2005), and particle velocities and entropies or internal energies are evolved according to a pressure-entropy or energy formalism (modern SPH) or a density-entropy or energy formalism (traditional SPH). Essentially, modern SPH techniques evaluate the pressure and the local density of each particle by considering the neighbours within the particle smoothing length and use pressure to define the equations of motion (Hopkins 2013). Traditional SPH techniques instead directly estimate the pressure of each particle from its local density in the same way as for the other particle properties, and use local density to define the equations of motion. In SPH methods the particle density mirrors the density of different regions of the fluid.

Grid-based and SPH techniques both have advantages, but also limits. At sufficiently high velocities, grid-based methods are non-invariant under Galilean transformations, which means that different reference frames are associated with different levels of numerical diffusion among adjacent cells, and simulation results may slightly depend on the choice of the reference system (e.g. Wadsley et al. 2008). Moreover, grid-based methods violate angular momentum conservation because a fluid moving across grid cells produces artificial diffusion; this diffusion can lead to unphysical forces, which couple with the fixed structure of the grid to tie the fluid motion on specific directions (e.g. Peery & Imlay 1988; Hahn et al. 2010). Finally, in grid-based methods hydrodynamics and gravity descriptions are mismatched, in the sense that hydrodynamics is evaluated by integrating quantities over each cell, while gravity is computed at the centre of each cell and then interpolated at the desired position (as for collisionless particles). This can produce spurious instabilities (e.g. Truelove et al. 1997).

SPH methods first need an artificial viscosity term added to the particle equation of motion in order to resolve shocks (Balsara 1989; Cullen & Dehnen 2010). Second, traditional SPH codes are associated with a surface tension between fluid regions of highly different densities, which limits their mixing (e.g. Agertz et al. 2007). Great effort has been made to improve SPH methods, leading to the so-called modern SPHs (Hopkins 2013). The smoothed definition of pressures together with densities, the more sophisticated viscosity switches, the higher order smoothing kernels (quintic spline instead of cubic spline; see below), and the inclusion of artificial conduction allowed solving these problems, at least partially. However, the higher order kernels typically lead to excessive diffusion. Despite all these improvements, some intrinsic limits of this technique still remain, such as the ideal infinite number of neighbours required to capture small-amplitude instabilities.

Recently, a completely new Lagrangian method that aims to simultaneously capture the advantages of both SPH and grid-based techniques, has been implemented in the public code GIZMO (Hopkins 2015). In GIZMO, the volume is discretised among a discrete set of tracers (particles) through a partition scheme based on a smoothing kernel (in a way that is similar to SPH codes). However, unlike SPH codes, these particles do

<footnote n="1">¹ http://www.mpa.mpa-garching.mpg.de/gadget/</footnote>

<footnote n="2">² http://www.tapir.caltech.edu/~phopkins/Site/GIZMO.html</footnote>

<footnote n="3">³ When a more realistic stellar equation of state is used (e.g. Rosswog et al. 2009, but only for white dwarfs), the value of $\beta_d$ may change slightly.</footnote>
<figure><img src="image_2.png" /><figcaption>**Fig. 1.** Left panel: analytic solutions for the $\gamma = 4/3, 5/3$ polytropic radial density profile from the Lane-Emden equation. Right panel: plot of the relaxed radial stellar density profile for each simulation technique for both polytropic indices ($\gamma = 4/3, 5/3$ from the highest to the lowest central density). Units are $M_\odot/R_\odot^2$ for $\rho$ and $R_\odot$ for $r$.</figcaption></figure>

not sample fluid elements, but only represent the centre of unstructured cells that are free to move with the fluid, like in moving mesh codes (Springel 2010). Hydrodynamics equations are then solved at the cell boundaries, defined by an effective face. This guarantees an exact conservation of energy and linear and angular momentum as well as an accurate description of shocks without needing an artificial viscosity term. The density associated with each particle or cell is obtained by dividing the mass of the cell for its effective volume. In this work, we use the mesh-free finite mass method of GIZMO, where particle mass is preserved, making the code perfectly Lagrangian. For this method, we use the cubic spline kernel with a desired number of neighbours equal to 32 for the partition.

### 3. SPH and GIZMO MFM simulations and stellar mass losses

We modelled stars as polytropes of index 5/3 (low-mass stars) or 4/3 (high-mass stars), with masses and radii of $1 M_\odot$ and $1 R_\odot$, sampling each of them with $N_{\text{part}} \sim 10^5$ particles. This is done by placing the particles through a close sphere packing and then stretching their radial positions to reach the required polytropic density profile, thus limiting the statistical noise associated with a random placement of the particles. $N_{\text{part}}$ sets the gravitational softening length of each particle in our codes to $\epsilon \sim 0.1 R_*/(N_{\text{part}})^{1/3} \sim 0.002 R_\odot$, preventing particle overlapping in evaluating gravitational interactions. We also tried test runs at higher resolution, where we modelled stars with $\sim 10^6$ particles, but did not find significant differences in the stellar mass loss $\Delta M$ with respect to simulations with lower resolution. We evolved stars in isolation for several dynamical times in order to ensure their stability. The right panel of Fig. 1 shows the relaxed stellar density profile, that is, the local density of the particles $\rho(r)$ (in $M_\odot/R_\odot^2$) versus their radial distance from the stellar centre of mass $r$ (in $R_\odot$), for each simulation technique for the two polytropic indices ($\gamma = 4/3$, and $5/3$ from the highest to the lowest central density), compared to the analytic solutions from the Lane-Emden equation (left panel). The kernel function that drives the evaluation of each particle local density (e.g. Price 2005) and the volume partition (Hopkins 2015) is chosen to be a cubic (in GADGET2 and GIZMO MFM) or quintic (in GIZMO modern SPH) spline, and the number of neighbours of each

particle and domain point within its smoothing length/kernel size is fixed to 32 and 128, respectively (Monaghan & Lattanzio 1985; Hongbin & Xin 2005; Dehnen & Aly 2012). Gravitational forces are computed through the Springel relative criterion (Springel 2005) instead of the standard Barnes-Hut criterion (Barnes & Hut 1986) because the Springel criterion shows better accuracy at the same computational cost. Since the relative criterion is based on the particle acceleration, which is not available at the beginning of each simulation, the Barnes-Hut criterion is adopted at the first timestep to estimate an acceleration value, and then the iteration is repeated using the Springel criterion in order to remain consistent with the subsequent iterations. In our simulations we use quite a large opening angle value (0.7), but the accuracy parameter for the relative criterion is set to 0.0005, which is very small compared to the suggested standard value (0.0025). We performed test runs setting the opening angle to 0.1 and increasing the accuracy parameter to 0.0025, but found no differences in the stellar density and temperature profiles and in $\Delta M$. We implemented the BH force through a Newtonian analytical potential, with $M_{\text{BH}} = 10^6 M_\odot$, and in each of the traditional SPH, modern SPH and GIZMO MFM simulations we placed one star on a parabolic orbit with a given $r_p$, that is, $\beta$, around the BH. The star was initially placed at a distance five times greater than $r_t$ to avoid spurious tidal distortions (we also tested larger initial distances, but found no significant differences in the outcomes). Stellar rotation is not expected to significantly affect our results in the range of $\beta$ considered in this paper (Stone et al. 2013). Figures 2 and 3 show snapshots from our traditional SPH simulations recorded shortly after pericentre passage. The lower limit of the range where $\beta_d$ lies (yellow) allows the core recollapse to occur for both polytropic indices (Guillochon & Ramirez-Ruiz 2013, 2015a). Modern SPH and GIZMO MFM simulations give almost the same results.

We aim to assess the stellar mass loss $\Delta M$ in each simulation. We recall that $\Delta M = M_*$ corresponds to total disruption. We describe the method we adopted to evaluate $\Delta M$ from each of our simulated star-BH tidal encounters, following GRR. In a specific simulation at a specific time, the position and velocity components of the stellar centre of mass around the BH, $x_{\text{CM}}$, $y_{\text{CM}}$, $z_{\text{CM}}$, $v_{x_{\text{CM}}}$, $v_{y_{\text{CM}}}$, and $v_{z_{\text{CM}}}$ are defined through an iterative approach. As a first step, we choose them to coincide with the position and velocity components of the particle with the highest
<figure><img src="image_3.png" /><figcaption>**Fig. 2.** Snapshots of the SPH particle density (in logarithmic scale) at $t \sim 8.5 \times 10^4$ s after pericentre passage for our GADGET2 simulations, in the case of a star with polytropic index 5/3. White and black correspond to the highest and lowest densities, respectively. Each snapshot is labelled with the corresponding value of $\beta$, with the range where the critical disruption parameter $\beta_d$ lies highlighted in yellow.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Fig. 3.** Same as Fig. 2 for a polytropic star of index 4/3.</figcaption></figure>
<figure><img src="image_5.png" /><figcaption>Fig. 4. Stellar mass loss (in units of ΔM/M*) as a function of β for a star with polytropic index 5/3. ΔM is evaluated at t ~ 10^6 s after the disruption. Blue, black, green, and red points are associated with GIZMO MFM, GADGET2, GIZMO modern SPH, and GRR simulations, respectively. Uncertainties on ΔM/M* from SPH and GIZMO MFM simulations are inferred as reported in the main text. Points at low values of β have been slightly horizontally displaced to give a better view of the error bars. The value of the critical disruption parameter β_d (dashed lines) slightly depends on the adopted simulation method.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>Fig. 5. Same as Fig. 4 for a polytropic star of index 4/3. The values of β_d obtained from our simulations visibly differ from those of GRR.</figcaption></figure>

local density, $x_{peak}$, $y_{peak}$, $z_{peak}$, $v_{x_{peak}}$, $v_{y_{peak}}$, $v_{z_{peak}}$. The specific binding energy to the star of the ith particle then reads

$$E_{*i} = \frac{1}{2}[(v_{x_i} - v_{x_{peak}})^2 + (v_{y_i} - v_{y_{peak}})^2 + (v_{z_i} - v_{z_{peak}})^2] + \phi_{*i}, \quad (4)$$

where $v_{x_i}$, $v_{y_i}$, and $v_{z_i}$ are the velocity components of the ith particle and $\phi_{*i}$ the stellar gravitational potential acting on the ith particle (directly computed by the simulation code). By considering only particles with $E_{*i} < 0$, we re-define the position and velocity components of the star centre of mass and re-evaluate Eq. (4) by setting them in place of the components labelled with

the subscript “peak”. The process is re-iterated until the convergence of $v_{CM}$ to a constant value to lower than $10^{-5} R_{\odot} \text{yr}^{-1}$. Particles with $E_{*i} > 0$ are unbound from the star. The stellar mass loss at the considered time can be obtained by multiplying the mass of a single particle, $m = M_{*}/N_{part}$, by the number of particles bound to the star, $N_{Bound}$, and subtracting the result from $M_{*}$. ΔM is obtained at t ~ 10^6 s (~650 stellar dynamical times) after the disruption.

Figures 4 and 5 show the stellar mass loss in units of ΔM/M* as a function of β for polytopes of index 5/3 and 4/3, respectively, inferred from our simulations with GIZMO MFM (blue
<table><caption>**Table 1.** Fitting coefficients of Eq. (6) and $\beta_d$ for each of the point sets in Figs. 4 and 5.</caption><thead><tr><th>Simulations</th><th>Polytropic index</th><th>A</th><th>B</th><th>C</th><th>D</th><th>E</th><th>β<sub>d</sub></th></tr></thead><tbody><tr><td>GRR</td><td>5/3</td><td>3.1647</td><td>-6.3777</td><td>3.1797</td><td>-3.4137</td><td>2.4616</td><td>0.90</td></tr><tr><td>GIZMO MFM</td><td>5/3</td><td>5.4722</td><td>-11.764</td><td>6.3204</td><td>-3.8172</td><td>2.8919</td><td>0.91</td></tr><tr><td>GADGET2</td><td>5/3</td><td>8.9696</td><td>-19.111</td><td>10.180</td><td>-4.2964</td><td>3.3231</td><td>0.93</td></tr><tr><td>GIZMO modern SPH</td><td>5/3</td><td>8.7074</td><td>-18.358</td><td>9.6760</td><td>-4.5340</td><td>3.5914</td><td>0.94</td></tr><tr><td>GRR</td><td>4/3</td><td>12.996</td><td>-31.149</td><td>12.865</td><td>-5.3232</td><td>6.4262</td><td>1.85</td></tr><tr><td>GIZMO MFM</td><td>4/3</td><td>-13.964</td><td>11.217</td><td>-2.1168</td><td>0.3930</td><td>0.5475</td><td>2.00</td></tr><tr><td>GADGET2</td><td>4/3</td><td>-15.378</td><td>-5.2385</td><td>6.3635</td><td>-1.5122</td><td>5.7378</td><td>2.02</td></tr><tr><td>GIZMO modern SPH</td><td>4/3</td><td>-10.394</td><td>-0.2160</td><td>2.6421</td><td>-0.8804</td><td>2.9215</td><td>2.02</td></tr></tbody></table>

<table><caption>**Table 2.** $\beta_d$ value as a function of polytropic index and adopted simulation method.</caption><thead><tr><th>Simulation method</th><th>Polytropic index</th><th>β<sub>d</sub></th></tr></thead><tbody><tr><td>AMR grid-based</td><td>5/3</td><td>0.90</td></tr><tr><td>MFM</td><td>5/3</td><td>0.91</td></tr><tr><td>Traditional SPH</td><td>5/3</td><td>0.93</td></tr><tr><td>Modern SPH</td><td>5/3</td><td>0.94</td></tr><tr><td>AMR grid-based</td><td>4/3</td><td>1.85</td></tr><tr><td>MFM</td><td>4/3</td><td>2.00</td></tr><tr><td>Traditional SPH</td><td>4/3</td><td>2.02</td></tr><tr><td>Modern SPH</td><td>4/3</td><td>2.02</td></tr></tbody></table>

points), GADGET2 (black points) and GIZMO modern SPH (green points), and the same obtained from the GRR simulations (red points). We estimate the uncertainty on our inferred $\Delta M/M_*$ as

$$ \sigma_{M*}^{\Delta M} = \sqrt{\sigma_{\text{Poisson}}^2 + \sigma_{E_{i}}^2 + \sigma_{\text{AD}}^2} = \sqrt{\left(\frac{\sqrt{N_{\text{Bound}}}}{N_{\text{part}}}\right)^2 + 0.01^2 + \sigma_{\text{AD}}^2} \quad (5) $$

where $\sigma_{\text{AD}}$ is the average deviation from 1 of $\Delta M/M_*$ for total disruptions in each of our point sets and $\sigma_{E_{i}} = 0.01$, as the values of $|E_{i}|$ for about $10^3$ particles of $10^5$ are lower than 0.01 times the average value $|E_{*}|$, that is, we are not able to determine exactly whether these $10^3$ particles are bound to or unbound from the star. We fit each of our point sets with a function introduced in GRR

$$ \begin{gathered} f(\beta) = \exp\left[\frac{A + B\beta + C\beta^2}{1 - D\beta + E\beta^2}\right], \quad \beta < \beta_d \\ f(\beta) = 1, \quad \beta \ge \beta_d. \end{gathered} \quad (6) $$

The values of the coefficients A-E and of $\beta_d$ are given in Table 1. It is worth noting that for the 5/3 polytropic index the curves of stellar mass loss associated with the four simulation codes differ very slightly in the value of the critical disruption parameter $\beta_d$ (dashed lines in Fig. 4). Specifically, $\beta_d$ is reached first in the GRR simulations ($\beta_d = 0.90$), followed by the GIZMO MFM ($\beta_d = 0.91$), GADGET2 ($\beta_d = 0.93$), and GIZMO modern SPH ($\beta_d = 0.94$) simulations (Table 2). This is expected because of the greater degree of excessive diffusion that characterises grid-based techniques compared to modern and

traditional SPH techniques and the surface tension conversely involved in SPH methods (Sect. 2). For the 4/3 polytropic index, instead, there is disagreement between our simulations and those of GRR (dashed lines in Fig. 5). $\beta_d$ is reached clearly first in the simulations of GRR ($\beta_d = 1.85$), followed by very similar values of the GIZMO MFM ($\beta_d = 2.00$), GADGET2 ($\beta_d = 2.02$), and GIZMO modern SPH ($\beta_d = 2.02$) simulations (Table 2). We hypothesise that the lower value of $\beta_d$ obtained by GRR is the result of resolving the stellar core of the $\gamma = 4/3$ polytrope not far enough.

In support of this hypothesis, we tested the dependence of $\beta_d$ on the resolution of our simulations by performing some low-resolution (~$10^3$ particles) GADGET2 simulations for the two polytropic indices (black points in Fig. 6). Figure 6 shows that for a $\gamma = 5/3$ polytrope (left-hand panel) the change in resolution has negligible effects on $\beta_d$. On the other hand, for $\gamma = 4/3$ polytopes (right-hand panel) we observe a strong dependence of $\beta_d$ on resolution below a resolution threshold because the configuration of the star is less stable.

We also determined the dependence of $\beta_d$ on different values of $M_{\text{BH}}$ by performing additional low-resolution (~$10^3$ particles) GADGET2 simulations with a $\gamma = 5/3$ polytrope of mass 1 $M_{\odot}$ and BHs of masses $10^5 M_{\odot}$ and $10^7 M_{\odot}$. Figure 7 clearly shows that $\beta_d$ does not depend sensitively on $M_{\text{BH}}$. We recall that flares and accretion temperatures instead depend on $M_{\text{BH}}$ (e.g. Guillochon & Ramirez-Ruiz 2013, 2015a).

For completeness, we also show in Fig. 8 how the polytropic index of the stellar remnant, which results from partial disruptions on parabolic orbits, is not preserved, but decreases with increasing $\beta$ for both $\gamma = 5/3$ polytopes (left panel) and $\gamma = 4/3$ polytopes (right panel).

## 4. Summary and conclusions

Tidal disruption events provide a unique way to probe otherwise quiescent or low-luminous black holes at the centres of galaxies. When approaching the central black hole of a galaxy, a star may be totally or partially disrupted by the black hole tidal field, depositing material onto the compact object and lighting it up through a bright flare (e.g. Rees 1988; Phinney 1989; Evans & Kochanek 1989). Such a tidal accretion flare is expected to be shaped by the structure of the disrupted star (e.g. Lodato et al. 2009) and the morphology of the star-black hole encounter (e.g. Guillochon & Ramirez-Ruiz 2013, 2015a).

The hydrodynamical simulations of Guillochon & Ramirez-Ruiz of star-black hole close encounters probably represent the most complete theoretical investigation of the
<figure><img src="image_7.png" /></figure>

**Fig. 6.** Comparison of mass losses as a function of $\beta$ near $\beta_d$ between the GRR simulations (blue points), high- (~10⁵ particles; red points) and low-resolution (~10³ particles; black points) GADGET2 simulations, for $\gamma = 5/3$ (left panel) and $\gamma = 4/3$ (right panel) polytropes. For a $\gamma = 4/3$ polytrope, the value of $\beta_d$ clearly depends on the adopted resolution below a resolution threshold. For a $\gamma = 5/3$ polytrope, the value of $\beta_d$ differs very slightly among the three simulations.

<figure><img src="image_8.png" /></figure>

**Fig. 7.** Comparison of mass losses as a function of $\beta$ near $\beta_d$ for a $\gamma = 5/3$ polytrope of mass $1 M_\odot$ approaching BHs with three different masses: $10^5 M_\odot$ (red points), $10^6 M_\odot$ (black points), $10^7 M_\odot$ (blue points). The value of $\beta_d$ clearly does not depend on $M_{BH}$.

properties of tidal disruption events (Guillochon & Ramirez-Ruiz 2013, 2015a). In each simulation, the star ($M_* = 1 M_\odot$, $R_* = 1 R_\odot$) is modelled as a polytrope of index 5/3 or 4/3 and evolved on a parabolic orbit with a specific pericentre around the black hole ($M_{BH} = 10^6 M_\odot$) using an AMR grid-based code. The resulting stellar mass loss defines the morphology of the simulated encounter, that is, it defines whether the disruption is total or partial, thus shaping the ensuing accretion flare. Here we followed the approach of Guillochon & Ramirez-Ruiz, but adopted two SPH simulation codes (GADGET2, traditional SPH; Springel 2005; GIZMO, modern SPH; Hopkins 2015) and GIZMO in MFM mode (Hopkins 2015) instead of a grid-based method, as all these simulation techniques have their advantages, but also limits (Sect. 2). We mainly intended to determine for each polytropic index whether the demarcation line between total and partial tidal disruption events, the critical disruption parameter $\beta_d$ (Eq. (3)), is the same for different simulation techniques.

Figures 4 and 5 clearly show that for a $\gamma = 5/3$ polytrope the curves of stellar mass loss inferred from AMR grid-based simulations (red points) and from GIZMO MFM (blue points), traditional SPH (black points), and modern SPH (green points) simulations differ only slightly in the value of $\beta_d$ (dashed lines), reflecting the limits of different codes (Sect. 2), while for a $\gamma = 4/3$ polytrope there is disagreement between our simulations and those of GRR (Table 2), which is most likely due to the adopted resolutions; this interpretation is consistent with the resolution tests we performed with our own simulations (Fig. 6). However, even with equal resolution, the SPH approach should be superior to a grid-based approach at resolving the dynamics of the core of, especially, a $\gamma = 4/3$ polytrope, given that the resolution naturally follows density in equal-mass-particle approaches. As a consequence, we find $\beta_d = 0.92 \pm 0.02$ ($2.01 \pm 0.01$) for a $\gamma = 5/3$ (4/3) polytrope.

The $\gamma = 4/3$ profile is probably only appropriate for a zero-age main-sequence sun because the central density of our Sun is about twice greater than the $\gamma = 4/3$ polytrope at an age of 5 Gyr. For a real star, even greater resolution would therefore be needed in a grid-based approach in order to properly estimate the location of full versus partial disruption. Moreover, real stars are generally not well modelled by a single polytropic index, especially as they evolve (MacLeod et al. 2012). Giant stars consist of a tenuous envelope and a dense core, which prevents envelope mass loss, thus likely moving the value of $\beta_d$ even ahead. A similar core-envelope structure and behaviour also characterise giant planets when they are disrupted by their host star (Liu et al. 2013). TDEs could also refer to disruptions by stellar objects (Guillochon et al. 2011; Perets et al. 2016). However, the value of $\beta_d$ for the latter encounters still remains to be investigated.

*Acknowledgements.* We thank the ISCRA staff for allowing us to perform our simulations on the Cineca Supercomputing Cluster GALILEO. We acknowledge P. F. Hopkins and G. Lodato for very useful discussion and comments on this work. E.R.C. acknowledges support by NASA through the Einstein Fellowship Program, grant PF6-170150. This work was also supported by the Packard grant and NASA ATP grant NNX14AH37G (ER). We also thank the referee, H. B. Perets, for valuable comments on the manuscript and very constructive suggestions.
<figure><img src="image_9.png" /><figcaption>**Fig. 8.** Changes in the value of the polytropic index of the stellar remnant resulting from partial disruptions for selected initial values of its *β*. Densities and radii are normalised to the central density and the radius of the remnant. Black curves represent solutions to the Lane-Emden equation for different values of *γ*: red, green, and blue points are from some of our simulations that left a remnant, for three different values of *β*. Left panel: *γ* = 5/3 polytrope. Right panel: *γ* = 4/3 polytrope.</figcaption></figure>

References

Agertz, O., Moore, B., Stadel, J., et al. 2007, MNRAS, 380, 963
Alexander, T. 2012, in EPJ Web of Conferences, 39, 5001
Balsara, D. S. 1989, Ph.D. Thesis, Univ. Illinois at Urbana-Champaign, USA
Barnes, J., & Hut, P. 1986, Nature, 324, 446
Berger, M. J., & Colella, P. 1989, J. Comput. Physics, 82, 64
Berger, M. J., & Oliger, J. 1984, J. Comput. Physics, 53, 484
Bonnerot, C., Rossi, E. M., Lodato, G., & Price, D. J. 2016, MNRAS, 455, 2253
Cannizzo, J. K., Lee, H. M., & Goodman, J. 1990, ApJ, 351, 38
Cossins, P. J. 2010, ArXiv e-prints [arXiv:1007.1245]
Coughlin, E. R., & Begelman, M. C. 2014, ApJ, 781, 82
Coughlin, E. R., & Nixon, C. 2015, ApJ, 808, 11
Coughlin, E. R., Nixon, C. et al. 2016, MNRAS, 455, 3612
Cullen, L., & Dehnen, W. 2010, MNRAS, 408, 669
Dehnen, W., & Aly, H. 2012, MNRAS, 425, 1068
Evans, C. R., & Kochanek, C. S. 1989, ApJ, 346, L13
Frank, J., & Rees, M. J. 1976, MNRAS, 176, 633
Guillochon, J., & Ramirez-Ruiz, E. 2013, ApJ, 767, 25
Guillochon, J., & Ramirez-Ruiz, E. 2015a, ApJ, 798, 64
Guillochon, J., & Ramirez-Ruiz, E. 2015b, ApJ, 809, 166
Guillochon, J., Ramirez-Ruiz, E., & Lin, D. 2011, ApJ, 732, 74
Hahn, O., Teyssier, R., & Carollo, M. 2010, MNRAS, 405, 274
Hayasaki, K., Stone, N., & Loeb, A. 2013, MNRAS, 434, 909
Hayasaki, K., Stone, N., & Loeb, A. 2016, MNRAS, 461, 3760
Hills, J. G. 1975, Nature, 254, 295
Ho, L. C. 2008, ARA&A, 46, 475
Hongbin, J., & Xin, D. 2005, J. Comput. Physics, 202, 699
Hopkins, P. F. 2013, MNRAS, 428, 2840
Hopkins, P. F. 2015, MNRAS, 450, 53
Kesden, M. 2012, Phys. Rev. D, 85, 024037
Kormendy, J., & Ho L. C. 2013, ARA&A, 51, 511
Kormendy, J., & Richstone, D. 1995, ARA&A, 33, 581
Lacy, J. H., Townes, C. H., & Hollenbach, D. J. 1982, ApJ, 262, 120
Laney, C. 1998, Computational Gasdynamics (Cambridge: Cambridge University Press)

Leveque, R. J. 1998, in Saas-Fee Advanced Course 27, Computational Methods for Astrophysical Fluid Flow (Berlin: Springer), eds. O. Steiner, & A. Gautschi, 1

Liu, S.-F., Guillochon, J. et al. 2013, ApJ, 762, 37

Lodato, G., & Rossi, E. M. 2011, MNRAS, 410, 359

Lodato, G., King, A. R., & Pringle, J. E. 2009, MNRAS, 392, 332

MacLeod, M., Guillochon, J., & Ramirez-Ruiz, E. 2012, ApJ, 757, 134

Magorrian, J., & Tremaine, S. 1999, MNRAS, 309, 447

Manukian, H., Guillochon, J. et al. 2013, ApJ, 771, L28

Monaghan, J. J. 1992, ARA&A, 30, 543

Monaghan, J. J., & Lattanzio, J. C. 1985, A&A, 149, 135

Peery, K., & Imlay, S. 1988, in Blunt-body flow simulations (American Institute of Aeronautics and Astronautics), 2904

Perets, H. B., Li, Z., Lombardi, J. C., Jr., & Milcarek, S. R., Jr. 2016, ApJ, 823, 113

Phinney, E. S. 1989, in The Center of the Galaxy, Proc. 136th IAU Symp., 543

Piran, T., Sadowsky, A., & Tchekhovskoy, A. 2015, MNRAS, 453, 157

Price, D. J. 2005, ArXiv e-print [arXiv:astro-ph/0507472]

Price, D. J. 2012, in Advances in Computational Astrophysics: Methods, Tools, and Outcome, eds. R. Capuzzo-Dolcetta, M. Limongi, & A. Tornambé, ASP Conf. Ser., 453, 249

Ramirez-Ruiz, E., & Rosswog, S. 2009, ApJ, 697, L77

Rees, M. J. 1988, Nature, 333, 523

Rosswog, S., Ramirez-Ruiz, E., & Hix, W. R. 2009, ApJ, 695, 404

Shiokawa, H., Krolik, J. H., Cheng, R. M., Piran, T., & Noble, S. C. 2015, ApJ, 804, 85

Springel, V. 2005, MNRAS, 364, 1105

Springel, V. 2010, MNRAS, 401, 791

Stone, N., Sari, R., & Loeb, A. 2013, MNRAS, 435, 1809

Strubbe, L. E., & Quataert, E. 2009, MNRAS, 400, 2070

Truelove, J. K., Klein, R. I. et al. 1997, ApJ, 489, L179

Vestergaard M., & Osmer P. S., 2009, ApJ, 699, 800

Wadsley, J. W., Veeravalli, G., & Couchman, H. M. P. 2008, MNRAS, 387, 427

Wang, J., & Merritt, D. 2004, ApJ, 600, 149