Issue 
A&A
Volume 600, April 2017



Article Number  A124  
Number of page(s)  8  
Section  Astrophysical processes  
DOI  https://doi.org/10.1051/00046361/201630092  
Published online  12 April 2017 
The fine line between total and partial tidal disruption events
^{1} Dipartimento di Fisica G. Occhialini, Università degli Studi di Milano Bicocca, Piazza della Scienza 3, 20126 Milano, Italy
email: d.mainetti1@campus.unimib.it
^{2} INAF, Osservatorio Astronomico di Brera, via E. Bianchi 46, 23807 Merate (LC), Italy
^{3} INFN, Sezione di MilanoBicocca, Piazza della Scienza 3, 20126 Milano, Italy
^{4} Institut d’Astrophysique de Paris, Sorbonne Universités, UPMC Univ Paris 6 et CNRS, UMR 7095, 98 bis bd Arago, 75014 Paris, France
^{5} Astronomy Department and Theoretical Astrophysics Center, University of California, Berkeley, CA 94720, USA
^{6} HarvardSmithsonian Center for Astrophysics, The Institute for Theory and Computation, 60 Garden Street, Cambridge, MA02138, USA
^{7} Department of Astronomy & Astrophysics, University of California, Santa Cruz, CA95064, USA
Received: 18 November 2016
Accepted: 24 February 2017
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 meshfree finite mass, traditional, and modern smoothed particle hydrodynamical simulations of starblack 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 socalled critical disruption parameter β_{d}) depends only weakly on the adopted simulation method, being β_{d} = 0.92 ± 0.02 for a γ = 5/3 polytrope and β_{d} = 2.01 ± 0.01 for a γ = 4/3 polytrope.
Key words: hydrodynamics / methods: numerical / galaxies: nuclei / black hole physics / accretion, accretion disks
© ESO, 2017
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 lowluminous 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 socalled BH tidal radius(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 & RamirezRuiz 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(2)where x encapsulates effects related to the BH spin (Kesden 2012). Hence, the nonrotating destroyer BH mass must be M_{BH} ≲ 10^{8}M_{⊙} when solartype stars are involved. TD accretion flares thus reveal otherwise quiescent or lowluminous 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 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; RamirezRuiz & Rosswog 2009; Hayasaki et al. 2013, 2016; Coughlin & Nixon 2015; Guillochon & RamirezRuiz 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 & RamirezRuiz 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 (3)for specific stellar structures that distinguishes total TDs from partial TDs, with r_{t} given by Eq. (1), β = r_{t}/r_{p} and r_{p} being the pericentre of the star around the BH. A partial TD is obtained for β<β_{d}, that is, for r_{p}>r_{d}, a total TD for β ≥ β_{d}, that is, for r_{p} ≤ 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 selfgravity only at the stellar surface, and not everywhere within the star. This problem has been considered previously (Guillochon & RamirezRuiz 2013, 2015a; hereafter GRR). GRR evaluated the critical disruption parameter β_{d} for polytropes of index 5/3 and 4/3 (which represent low and highmass stars, respectively) using a series of adaptive mesh refinement (AMR) gridbased 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)^{1} and gizmo (modern SPH and meshfree finite mass (MFM); Hopkins 2015)^{2}. Since these techniques all have advantages but also limits, we are inclined to follow GRR in finding the critical disruption parameter β_{d} for certain stellar structures^{3} using an MFM, a traditional SPH, and a modern SPH code instead of an AMR gridbased code, with the goal of comparing results from different techniques.
The paper is organised as follows. In Sect. 2 we compare AMR gridbased 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 ΔM in our simulated encounters. We show the curves of mass loss we obtain for all codes as a function of β and polytropic index, comparing them and the corresponding β_{d} with GRR. Section 4 contains our summary.
2. Gridbased vs. SPH and gizmo mfm codes
Fluid hydrodynamics and interactions in astrophysics are generally treated using two different classes of numerical methods: Eulerian gridbased (e.g. Laney 1998; Leveque 1998) and Lagrangian SPH (e.g. Monaghan 1992; Price 2005, 2012; Cossins 2010). Basically, gridbased 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 gridbased 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 highdensity 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 socalled smoothing length (e.g. Price 2005), and particle velocities and entropies or internal energies are evolved according to a pressureentropy or energy formalism (modern SPH) or a densityentropy 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.
Gridbased and SPH techniques both have advantages, but also limits. At sufficiently high velocities, gridbased methods are noninvariant 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, gridbased 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 gridbased 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 socalled 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 smallamplitude instabilities.
Fig. 1 Left panel: analytic solutions for the γ = 4/3,5/3 polytropic radial density profile from the LaneEmden equation. Right panel: plot of the relaxed radial stellar density profile for each simulation technique for both politropic indices (γ = 4/3,5/3 from the highest to the lowest central density). Units are for ρ and R_{⊙} for r. 
Recently, a completely new Lagrangian method that aims to simultaneously capture the advantages of both SPH and gridbased 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 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 meshfree 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
Fig. 2 Snapshots of the SPH particle density (in logarithmic scale) at t ~ 8.5 × 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 β, with the range where the critical disruption parameter β_{d} lies highlighted in yellow. 
We modelled stars as polytropes of index 5/3 (lowmass stars) or 4/3 (highmass stars), with masses and radii of 1 M_{⊙} and 1 R_{⊙}, sampling each of them with N_{part} ~ 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_{part} sets the gravitational softening length of each particle in our codes to ϵ ~ 0.1 R_{∗}/ (N_{part})^{1/3} ~ 0.002 R_{⊙}, preventing particle overlapping in evaluating gravitational interactions. We also tried test runs at higher resolution, where we modelled stars with ~10^{6} particles, but did not find significant differences in the stellar mass loss Δ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 ρ(r) (in ) versus their radial distance from the stellar centre of mass r (in R_{⊙}), for each simulation technique for the two polytropic indices (γ = 4/3, and 5/3 from the highest to the lowest central density), compared to the analytic solutions from the LaneEmden 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 BarnesHut 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 BarnesHut 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 ΔM. We implemented the BH force through a Newtonian analytical potential, with M_{BH} = 10^{6}M_{⊙}, 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, β, 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 β 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 β_{d} lies (yellow) allows the core recollapse to occur for both polytropic indices (Guillochon & RamirezRuiz 2013, 2015a). Modern SPH and gizmo mfm simulations give almost the same results.
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. 
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. 
We aim to assess the stellar mass loss ΔM in each simulation. We recall that ΔM = M_{∗} corresponds to total disruption. We describe the method we adopted to evaluate ΔM from each of our simulated starBH 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_{CM}, y_{CM}, z_{CM}, v_{xCM}, v_{yCM}, and v_{zCM} 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 local density, x_{peak}, y_{peak}, z_{peak}, v_{xpeak}, v_{ypeak}, v_{zpeak}. The specific binding energy to the star of the ith particle then reads (4)where v_{xi}, v_{yi}, and v_{zi} are the velocity components of the ith particle and φ_{∗ 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 redefine the position and velocity components of the star centre of mass and reevaluate Eq. (4) by setting them in place of the components labelled with the subscript “peak”. The process is reiterated until the convergency of v_{CM} to a constant value to lower than 10^{5}R_{⊙} 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.
β_{d} value as a function of polytropic index and adopted simulation method.
Figures 4 and 5 show the stellar mass loss in units of ΔM/M_{∗} as a function of β for polytropes of index 5/3 and 4/3, respectively, inferred from our simulations with gizmo mfm (blue 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 ΔM/M_{∗} as (5)where σ_{AD} is the average deviation from 1 of ΔM/M_{∗} for total disruptions in each of our point sets and σ_{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 , 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 (6)The values of the coefficients A−E and of β_{d} are given in Table 1.
Fig. 6 Comparison of mass losses as a function of β near β_{d} between the GRR simulations (blue points), high (~10^{5} particles; red points) and lowresolution (~10^{3} particles; black points) gadget2 simulations, for γ = 5/3 (left panel) and γ = 4/3 (right panel) polytropes. For a γ = 4/3 polytrope, the value of β_{d} clearly depends on the adopted resolution below a resolution threshold. For a γ = 5/3 polytrope, the value of β_{d} differs very slightly among the three simulations. 
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 β_{d} (dashed lines in Fig. 4). Specifically, β_{d} is reached first in the GRR simulations (β_{d} = 0.90), followed by the gizmo mfm (β_{d} = 0.91), gadget2 (β_{d} = 0.93), and gizmo modern SPH (β_{d} = 0.94) simulations (Table 2). This is expected because of the greater degree of excessive diffusion that characterises gridbased 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). β_{d} is reached clearly first in the simulations of GRR (β_{d} = 1.85), followed by very similar values of the gizmo mfm (β_{d} = 2.00), gadget2 (β_{d} = 2.02), and gizmo modern SPH (β_{d} = 2.02) simulations (Table 2). We hypothesise that the lower value of β_{d} obtained by GRR is the result of resolving the stellar core of the γ = 4/3 polytrope not far enough.
In support of this hypothesis, we tested the dependence of β_{d} on the resolution of our simulations by performing some lowresolution (~10^{3} particles) gadget2 simulations for the two polytropic indices (black points in Fig. 6). Figure 6 shows that for a γ = 5/3 polytrope (lefthand panel) the change in resolution has negligible effects on β_{d}. On the other hand, for γ = 4/3 polytropes (righthand panel) we observe a strong dependence of β_{d} on resolution below a resolution threshold because the configuration of the star is less stable.
We also determined the dependence of β_{d} on different values of M_{BH} by performing additional lowresolution (~10^{3} particles) gadget2 simulations with a γ = 5/3 polytrope of mass 1 M_{⊙} and BHs of masses 10^{5}M_{⊙} and 10^{7}M_{⊙}. Figure 7 clearly shows that β_{d} does not depend sensitively on M_{BH}. We recall that flares and accretion temperatures instead depend on M_{BH} (e.g. Guillochon & RamirezRuiz 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 β for both γ = 5/3 polytropes (left panel) and γ = 4/3 polytropes (right panel).
Fig. 7 Comparison of mass losses as a function of β near β_{d} for a γ = 5/3 polytrope of mass 1 M_{⊙} approaching BHs with three different masses: 10^{5}M_{⊙} (red points), 10^{6}M_{⊙} (black points), 10^{7}M_{⊙} (blue points). The value of β_{d} clearly does not depend on M_{BH}. 
4. Summary and conclusions
Tidal disruption events provide a unique way to probe otherwise quiescent or lowluminous 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 starblack hole encounter (e.g. Guillochon & RamirezRuiz 2013, 2015a).
The hydrodynamical simulations of Guillochon & RamirezRuiz of starblack hole close encounters probably represent the most complete theoretical investigation of the properties of tidal disruption events (Guillochon & RamirezRuiz 2013, 2015a). In each simulation, the star (M_{∗} = 1 M_{⊙}, R_{∗} = 1 R_{⊙}) 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_{⊙}) using an AMR gridbased 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 & RamirezRuiz, 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 gridbased 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 β_{d} (Eq. (3)), is the same for different simulation techniques.
Figures 4 and 5 clearly show that for a γ = 5/3 polytrope the curves of stellar mass loss inferred from AMR gridbased 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 β_{d} (dashed lines), reflecting the limits of different codes (Sect. 2), while for a γ = 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 gridbased approach at resolving the dynamics of the core of, especially, a γ = 4/3 polytrope, given that the resolution naturally follows density in equalmassparticle approaches. As a consequence, we find β_{d} = 0.92 ± 0.02 (2.01 ± 0.01) for a γ = 5/3 (4/3) polytrope.
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 LaneEmden 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. 
The γ = 4/3 profile is probably only appropriate for a zeroage mainsequence sun because the central density of our Sun is about twice greater than the γ = 4/3 polytrope at an age of 5 Gyr. For a real star, even greater resolution would therefore be needed in a gridbased 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 β_{d} even ahead. A similar coreenvelope 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 β_{d} for the latter encounters still remains to be investigated.
When a more realistic stellar equation of state is used (e.g. Rosswog et al. 2009, but only for white dwarfs), the value of β_{d} may change slightly.
Acknowledgments
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 PF6170150. 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.
References
 Agertz, O., Moore, B., Stadel, J., et al. 2007, MNRAS, 380, 963 [NASA ADS] [CrossRef] [Google Scholar]
 Alexander, T. 2012, in EPJ Web of Conferences, 39, 5001 [Google Scholar]
 Balsara, D. S. 1989, Ph.D. Thesis, Univ. Illinois at UrbanaChampaign, USA [Google Scholar]
 Barnes, J., & Hut, P. 1986, Nature, 324, 446 [NASA ADS] [CrossRef] [Google Scholar]
 Berger, M. J., & Colella, P. 1989, J. Comput. Physics, 82, 64 [NASA ADS] [CrossRef] [Google Scholar]
 Berger, M. J., & Oliger, J. 1984, J. Comput. Physics, 53, 484 [Google Scholar]
 Bonnerot, C., Rossi, E. M., Lodato, G., & Price, D. J. 2016, MNRAS, 455, 2253 [NASA ADS] [CrossRef] [Google Scholar]
 Cannizzo, J. K., Lee, H. M., & Goodman, J. 1990, ApJ, 351, 38 [NASA ADS] [CrossRef] [Google Scholar]
 Cossins, P. J. 2010, ArXiv eprints [arXiv:1007.1245] [Google Scholar]
 Coughlin, E. R., & Begelman, M. C. 2014, ApJ, 781, 82 [NASA ADS] [CrossRef] [Google Scholar]
 Coughlin, E. R., & Nixon, C. 2015, ApJ, 808, 11 [NASA ADS] [CrossRef] [Google Scholar]
 Coughlin, E. R., Nixon, C., et al. 2016, MNRAS, 455, 3612 [NASA ADS] [CrossRef] [Google Scholar]
 Cullen, L., & Dehnen, W. 2010, MNRAS, 408, 669 [NASA ADS] [CrossRef] [Google Scholar]
 Dehnen, W., & Aly, H. 2012, MNRAS, 425, 1068 [NASA ADS] [CrossRef] [Google Scholar]
 Evans, C. R., & Kochanek, C. S. 1989, ApJ, 346, L13 [NASA ADS] [CrossRef] [Google Scholar]
 Frank, J., & Rees, M. J. 1976, MNRAS, 176, 633 [NASA ADS] [CrossRef] [Google Scholar]
 Guillochon, J., & RamirezRuiz, E. 2013, ApJ, 767, 25 [NASA ADS] [CrossRef] [Google Scholar]
 Guillochon, J., & RamirezRuiz, E. 2015a, ApJ, 798, 64 [NASA ADS] [CrossRef] [Google Scholar]
 Guillochon, J., & RamirezRuiz, E. 2015b, ApJ, 809, 166 [NASA ADS] [CrossRef] [Google Scholar]
 Guillochon, J., RamirezRuiz, E., & Lin, D. 2011, ApJ, 732, 74 [NASA ADS] [CrossRef] [Google Scholar]
 Hahn, O., Teyssier, R., & Carollo, M. 2010, MNRAS, 405, 274 [NASA ADS] [Google Scholar]
 Hayasaki, K., Stone, N., & Loeb, A. 2013, MNRAS, 434, 909 [NASA ADS] [CrossRef] [Google Scholar]
 Hayasaki, K., Stone, N., & Loeb, A. 2016, MNRAS, 461, 3760 [NASA ADS] [CrossRef] [Google Scholar]
 Hills, J. G. 1975, Nature, 254, 295 [NASA ADS] [CrossRef] [Google Scholar]
 Ho, L. C. 2008, ARA&A, 46, 475 [NASA ADS] [CrossRef] [Google Scholar]
 Hongbin, J., & Xin, D. 2005, J. Comput. Physics, 202, 699 [NASA ADS] [CrossRef] [Google Scholar]
 Hopkins, P. F. 2013, MNRAS, 428, 2840 [NASA ADS] [CrossRef] [Google Scholar]
 Hopkins, P. F. 2015, MNRAS, 450, 53 [NASA ADS] [CrossRef] [Google Scholar]
 Kesden, M. 2012, Phys. Rev. D, 85, 024037 [NASA ADS] [CrossRef] [Google Scholar]
 Kormendy, J., & Ho, L. C. 2013, ARA&A, 51, 511 [Google Scholar]
 Kormendy, J., & Richstone, D. 1995, ARA&A, 33, 581 [NASA ADS] [CrossRef] [Google Scholar]
 Lacy, J. H., Townes, C. H., & Hollenbach, D. J. 1982, ApJ, 262, 120 [NASA ADS] [CrossRef] [Google Scholar]
 Laney, C. 1998, Computational Gasdynamics (Cambridge: Cambridge University Press) [Google Scholar]
 Leveque, R. J. 1998, in SaasFee Advanced Course 27, Computational Methods for Astrophysical Fluid Flow (Berlin: Springer), eds. O. Steiner, & A. Gautschy, 1 [Google Scholar]
 Liu, S.F., Guillochon, J., et al. 2013, ApJ, 762, 37 [NASA ADS] [CrossRef] [Google Scholar]
 Lodato, G., & Rossi, E. M. 2011, MNRAS, 410, 359 [NASA ADS] [CrossRef] [Google Scholar]
 Lodato, G., King, A. R., & Pringle, J. E. 2009, MNRAS, 392, 332 [NASA ADS] [CrossRef] [Google Scholar]
 MacLeod, M., Guillochon, J., & RamirezRuiz, E. 2012, ApJ, 757, 134 [NASA ADS] [CrossRef] [Google Scholar]
 Magorrian, J., & Tremaine, S. 1999, MNRAS, 309, 447 [NASA ADS] [CrossRef] [Google Scholar]
 Manukian, H., Guillochon, J., et al. 2013, ApJ, 771, L28 [NASA ADS] [CrossRef] [Google Scholar]
 Monaghan, J. J. 1992, ARA&A, 30, 543 [NASA ADS] [CrossRef] [Google Scholar]
 Monaghan, J. J., & Lattanzio, J. C. 1985, A&A, 149, 135 [NASA ADS] [Google Scholar]
 Peery, K., & Imlay, S. 1988, in Bluntbody flow simulations (American Institute of Aeronautics and Astronautics), 2904 [Google Scholar]
 Perets, H. B., Li, Z., Lombardi, J. C., Jr., & Milcarek, S. R., Jr. 2016, ApJ, 823, 113 [NASA ADS] [CrossRef] [Google Scholar]
 Phinney, E. S. 1989, in The Center of the Galaxy, Proc. 136th IAU Symp., 543 [Google Scholar]
 Piran, T., Sadowsky, A., & Tchekhovskoy, A. 2015, MNRAS, 453, 157 [NASA ADS] [CrossRef] [Google Scholar]
 Price, D. J. 2005, ArXiv eprint [arXiv:astroph/0507472] [Google Scholar]
 Price, D. J. 2012, in Advances in Computational Astrophysics: Methods, Tools, and Outcome, eds. R. CapuzzoDolcetta, M. Limongi, & A. Tornambé, ASP Conf. Ser., 453, 249 [Google Scholar]
 RamirezRuiz, E., & Rosswog, S. 2009, ApJ, 697, L77 [NASA ADS] [CrossRef] [Google Scholar]
 Rees, M. J. 1988, Nature, 333, 523 [NASA ADS] [CrossRef] [Google Scholar]
 Rosswog, S., RamirezRuiz, E., & Hix, W. R. 2009, ApJ, 695, 404 [NASA ADS] [CrossRef] [Google Scholar]
 Shiokawa, H., Krolik, J. H., Cheng, R. M., Piran, T., & Noble, S. C. 2015, ApJ, 804, 85 [NASA ADS] [CrossRef] [Google Scholar]
 Springel, V. 2005, MNRAS, 364, 1105 [NASA ADS] [CrossRef] [Google Scholar]
 Springel, V. 2010, MNRAS, 401, 791 [NASA ADS] [CrossRef] [Google Scholar]
 Stone, N., Sari, R., & Loeb, A. 2013, MNRAS, 435, 1809 [NASA ADS] [CrossRef] [Google Scholar]
 Strubbe, L. E., & Quataert, E. 2009, MNRAS, 400, 2070 [NASA ADS] [CrossRef] [Google Scholar]
 Truelove, J. K., Klein, R. I., et al. 1997, ApJ, 489, L179 [NASA ADS] [CrossRef] [Google Scholar]
 Vestergaard, M., & Osmer, P. S., 2009, ApJ, 699, 800 [Google Scholar]
 Wadsley, J. W., Veeravalli, G., & Couchman, H. M. P. 2008, MNRAS, 387, 427 [NASA ADS] [CrossRef] [Google Scholar]
 Wang, J., & Merritt, D. 2004, ApJ, 600, 149 [NASA ADS] [CrossRef] [Google Scholar]
All Tables
All Figures
Fig. 1 Left panel: analytic solutions for the γ = 4/3,5/3 polytropic radial density profile from the LaneEmden equation. Right panel: plot of the relaxed radial stellar density profile for each simulation technique for both politropic indices (γ = 4/3,5/3 from the highest to the lowest central density). Units are for ρ and R_{⊙} for r. 

In the text 
Fig. 2 Snapshots of the SPH particle density (in logarithmic scale) at t ~ 8.5 × 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 β, with the range where the critical disruption parameter β_{d} lies highlighted in yellow. 

In the text 
Fig. 3 Same as Fig. 2 for a polytropic star of index 4/3. 

In the text 
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. 

In the text 
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. 

In the text 
Fig. 6 Comparison of mass losses as a function of β near β_{d} between the GRR simulations (blue points), high (~10^{5} particles; red points) and lowresolution (~10^{3} particles; black points) gadget2 simulations, for γ = 5/3 (left panel) and γ = 4/3 (right panel) polytropes. For a γ = 4/3 polytrope, the value of β_{d} clearly depends on the adopted resolution below a resolution threshold. For a γ = 5/3 polytrope, the value of β_{d} differs very slightly among the three simulations. 

In the text 
Fig. 7 Comparison of mass losses as a function of β near β_{d} for a γ = 5/3 polytrope of mass 1 M_{⊙} approaching BHs with three different masses: 10^{5}M_{⊙} (red points), 10^{6}M_{⊙} (black points), 10^{7}M_{⊙} (blue points). The value of β_{d} clearly does not depend on M_{BH}. 

In the text 
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 LaneEmden 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. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.