In-situ extreme mass ratio inspirals via sub-parsec formation and migration of stars in thin, gravitationally unstable AGN discs

We investigate the properties of stars born via gravitational instability in accretion discs around supermassive black holes (SMBHs) in active galactic nuclei (AGN), and how this varies with the SMBH mass, accretion rate, or viscosity. We show with geometrically thin, steady-state disc solutions that fragmentation results in different populations of stars when one considers the initial conditions (e.g. density and temperature of the gravitationally unstable regions). We find that opacity gaps in discs around $10^6 \rm M_{\odot}$ SMBHs can trigger fragmentation at radii $\lesssim 10^{-2}$ pc, although the conditions lead to the formation of initially low stellar masses around $0.1-0.5 \rm M_{\odot}$. Discs around more massive SMBHs ($M_{\rm BH} =10^{7-8} \rm M_{\odot}$) form moderately massive or supermassive stars (the majority at $10^{0-2} \rm M_{\odot}$). Using linear migration estimates, we discuss three outcomes: stars migrate till they are tidally destroyed, accreted as extreme mass ratio inspirals (EMRIs), or leftover after disc dispersal. For a single AGN activity cycle, we find a lower-limit for the EMRI rate $R_{\rm emri}\sim 0-10^{-4} \rm yr^{-1}$ per AGN assuming a SF efficiency $\epsilon=1-30\%$. In cases where EMRIs occur, this implies a volumetric rate up to $0.5-10 \rm yr^{-1} Gpc^{-3}$ in the local Universe. The rates are particularly sensitive to model parameters for $M_{\rm BH}=10^6 \rm M_{\odot}$, for which EMRIs only occur if stars can accrete to $10$s of solar masses. Our results provide further evidence that gas-embedded EMRIs can contribute a substantial fraction of events detectable by milliHz gravitational wave detectors such as LISA. Our disc solutions suggest the presence of migration traps, as has been found for more massive SMBH discs. Finally, the surviving population of stars after the disc lifetime leaves implications for stellar discs in galactic nuclei.


INTRODUCTION
The most popular interpretation of quasars proposes that supermassive black holes (SMBHs) with masses BH 10 5 M reside in galactic nuclei surrounded by dense gaseous accretion discs. The bright emission originates from the accretion of gas onto the SMBH, which is driven by some mechanism of viscous dissipation. While the precise structure of the accretion flow onto SMBHs is poorly constrained and may take many forms, models of geometrically thin, near-Keplerian discs (e.g. Shakura & Sunyaev 1973) reproduce observed continuum emission in several cases leading such geometry to be considered a standard benchmark for accretion flow in quasars (Krolik 1999). Similar models explain active galactic nuclei (AGN) throughout the Universe to varying degrees, particularly those that are bright or shining at near-Eddington rates. An inevitable implication of thin disc models, however, is their vulnerability to gravitational instability in the outer regions (e.g. Shlosman & Begelman 1987;Goodman 2003).
Gravitational instability (GI), initially introduced in the context of stellar discs by Toomre (1964), occurs when the self-gravity of the gas overcomes its pressure and rotational support. Depending on the ability of the gas to cool, the gravitationally unstable disc can fragment into coherent structures, leading to the formation of bound clumps (Gammie 2001), possibly halting the accretion flow entirely ★ E-mail: andrea@ics.uzh.ch Shlosman et al. 1990;Inayoshi & Haiman 2016). GI also generates spiral density waves (Durisen et al. 2007). These nonaxisymmetries can dominate the angular momentum transport in heavily self-gravitating discs (e.g. Kumar 1999;Cossins et al. 2009;Lodato & Rice 2004, 2005 and ultimately determine the disc evolution (Boley et al. 2010). In turn, energy dissipation due to gas collapse and density wave propagation generates heat, which can re-stabilise the disc. The interplay between collapse and restabilization suggests that regions of discs may exist in a quasisteady state of inflow, or a gravitoturbulent state in the case that fragmentation is suppressed (Rafikov 2009), where the accretion flow is partially driven albeit not terminated by GI.
Fragmentation due to marginal gravitational instability has been proposed as a mechanism for star formation in AGN discs (Shlosman & Begelman 1987;Wang & Silk 1994;Levin 2007;Goodman 2003;Goodman & Tan 2004;Nayakshin 2006) that may explain stellar discs in the Galactic center (Levin & Beloborodov 2003;Nayakshin & Sunyaev 2005;Davies & Lin 2020), in addition to altering the growth rate of SMBHs (Dittmann & Miller 2020). Analogously it is considered to be a formation mechanism for planets in gravitationally unstable protoplanetary discs (Boss 1997;Helled et al. 2014). Star formation in AGN discs and subsequent growth and interactions of embedded compact remnants (CRs) can lead to interesting observational signatures in both electromagnetic emission and gravitational waves (GWs) across a range of frequencies. AGN-assisted CR in-teraction is a proposed formation channel for high-frequency GW events detectable by LIGO-Virgo, in particular to explain detected stellar origin BH-BH mergers with large constituent masses and events with unequal mass ratios or unusual spins (Bartos et al. 2017;Stone et al. 2017;Secunda et al. 2019;Yang et al. 2019;McKernan et al. 2018McKernan et al. , 2020Tagawa et al. 2020aTagawa et al. ,b, 2021Gayathri et al. 2021). Embedded CRs (or stars, if not disrupted) may also migrate to the central regions of the disc where they interact with the central SMBH and become lower-frequency GW sources, referred to as extreme-or intermediate-mass ratio inspirals (E/IMRIs), which will be detectable by future space-based interferometers such as LISA (Amaro-Seoane et al. 2017), TianQin (Fan et al. 2020), or TAĲI Gong et al. (2021).
EMRIs are a prime target for milliHz GW missions. The detection of even one will provide high precision measurements of the central SMBH mass and spin (Barack & Cutler 2007;Babak et al. 2017). If the EMRI is embedded in gas, the GW signal can carry signatures of the environment that can be potentially constrain environment properties and complex migration physics Kocsis et al. 2011;Barausse et al. 2014;Derdzinski et al. 2019Derdzinski et al. , 2021Zwick et al. 2021). Rate estimates of EMRI formation via dynamical interactions in galactic nuclei span orders of magnitude in uncertainty (Amaro-Seoane 2020). Luckily, it is becoming apparent that the presence of an accretion disc can boost the EMRI formation rate . One mechanism for EMRI facilitation occurs due to the influence of the accretion disc on the orbits of stars and CRs in the nucleus, which tend to align and corotate with the disc (Syer et al. 1991;Norman & Silk 1983;Fabj et al. 2020;Tagawa et al. 2020a). Additionally, instability-induced star formation can produce in-situ CRs that accrete, merge, or migrate towards the central SMBH (Levin 2007;Goodman & Tan 2004). The latter works describe the collapse and aggregation of clumps in marginally stable discs that leads to the formation of a single, supermassive embedded star.
In this work, we revisit models of gravitationally unstable accretion discs, extending solutions over a set of system parameters (namely SMBH mass and rate of viscous inflow), to demonstrate the range of initial conditions for in-situ AGN stars. We investigate the diversity of the stellar mass distribution and subsequent evolutionary outcomes. Not all AGN are the same-SMBH masses, the total gas mass, and efficiency of inflow can vary, leading to different predictions for the location of fragmentation and the AGN star evolutionary outcomes: tidal disruption, in-situ EMRIs, or leftover stellar populations. In the following Section, we review recent developments of fragmentation in the protoplanetary disc context, for which several numerical investigations find deviations from the standard theory when attempting to explain properties of planets formed through GI. We then briefly summarize theoretical aspects of SF employed in this work, before describing our disc model in Section 2.2. With solutions for disc structure, we derive the regime where fragmentation is expected to occur in Section 2.3. In Section 3 we describe the process of fragmentation to gravitational collapse of protostars, based on the size of the instability, the collapse rate, and the accretion onto the resulting cores, and we derive the initial mass distribution of stars. In Section 4 we describe the subsequent evolution of AGN-star populations, given their ability to migrate and evolve throughout the disc. In Section 5 we discuss the outcomes of embedded stars, which can be tidally disrupted or accreted by the SMBH, and we quantify the rate of in-situ EMRIs for each disc solution. Finally, in Section 6 we summarise the takeaways and caveats of this work.

Astrophysical rationale
The criteria for fragmentation of gas into bound objects has been extensively tested in numerical simulations, often in the context of protoplanetary discs. It is well understood that fragmentation will occur if the disc cools on a timescale comparable to the local orbital time (Gammie 2001;Mayer & Gawryszczak 2008;Johnson & Gammie 2003;Clarke et al. 2007), although determining the precise conditions has proved a challenging numerical problem. The criteria for fragmentation, as well as the outcome of collapsing clumps (whether at the scale of planets, stars, or something in between) are sensitive to numerical implementations of disc thermodynamics Mayer et al. 2004), numerical viscosity (Hu et al. 2014), gravitational softening (Young & Clarke 2015), and small-scale resolution (Nelson 2006;Paardekooper 2012;Brucy & Hennebelle 2021). Recent numerical studies find convergence (over implementations of artificial viscosity) on the critical cooling rate below which fragmentation occurs, finding crit ≡ cool Ω ≈ 3 (Deng et al. 2017), where cool is the radiative cooling timescale and Ω = ( / 3 ) 1/2 is the orbital frequency at a radius in a disc around a central mass . This is the criteria we adopt in this work.
Convergence on the size and mass distribution of resulting clumps has yet to be attained numerically, however. Linear perturbation theory that computes the growth rate of local axisymmetric density perturbations in razor-thin discs predicts that fragmentation produces clumps on a characteristic scale, namely the most unstable Toomre wavelength (Toomre 1964). Yet several simulations find that clumps form at systematically smaller sizes, even differing by an order of magnitude from the predictions of perturbation theory (e.g. Boley et al. 2010;Galvagni et al. 2012;Galvagni & Mayer 2014;Müller et al. 2018;Deng et al. 2021;Forgan & Rice 2011). Differences arise due to finite disc thickness, departure from axisymmetry, and nonlinear effects in the gas flow dynamics immediately prior to fragmentation. In particular, fragmentation is almost invariably seen to occur along spiral arms generated in the initial phase of a gravitational instability in realistic discs on a variety of scales, hence clumps form in a gas flow having inherently a non-axisymmetric character (see, e.g. Boley et al. 2010 for protoplanetary disc scales or Tamburello et al. 2015 for galaxy-scale fragmentation). In addition, soon after fragmentation ensues, other effects matter in the nonlinear stage, which, by construction, cannot be captured by perturbation theory, such as the rotational support in clumps (Galvagni et al. 2012), and their vulnerability by local shear and stellar tidal forces (Galvagni & Mayer 2014;Müller et al. 2018). Furthermore, and perhaps most critical for ionised gas in AGN discs, large deviations from predictions of Toomre instability theory, both before and after fragmentation, arise in simulations of discs threaded with magnetic fields. This results in clump sizes that are even two orders of magnitude smaller than those in non-magnetized discs (Deng et al. 2020). When estimating the characteristic sizes of pre-stellar clumps in this work, we draw from the results of recent numerical simulations to take these reductions into account.
Considering that we are concerned with fragmentation at the AGN disc scale, we must also consider basic details of star formation. Conventional star formation involves the collapse of diffuse gas to dense, accreting cores, involving a vast range of temperatures and physical scales. Determining the stellar mass distribution from a collapsing cloud is a non-trivial problem given that complexities in opacity, radiative processes, rotation, and magnetic fields are all critical for determining the resulting stellar masses (Larson 1969;Penston 1969;Stahler & Palla 2004), particularly in regimes of high accretion rates that lead to supermassive stars (Haemmerlé et al. 2021;Hosokawa 2019).
In the case of star formation induced by GI in an AGN disc, fragmentation occurs in high angular momentum material and at initially higher densities and temperatures than a typical molecular cloud. Protostars effectively skip the initially slow phase of prestellar contraction 1 . Such cores collapse rapidly due to the gas conditions, as we will show, experiencing high accretion rates. Certain conditions lead to the formation of supermassive stars (SMSs), for which there are greater uncertainties in the final stellar masses and structure (Hosokawa et al. 2013) and subsequent evolution (Haemmerlé et al. 2021). Critically, the formation and subsequent evolution of AGNstars does not occur in isolation: gas from the disc may feed the protostellar envelopes and exert torques.

Disc model
We define a system of equations describing a thin, Keplerian, steadystate accretion disc. The central BH accretes at a constant fraction of the Eddington rate = Edd Edd , where the Eddington rate Edd = 4 BH eff es for a central BH of mass BH and a radiative efficiency eff = 0.1. The viscosity is parameterised with the alphadisc prescription, as in Shakura & Sunyaev (1973). The disc model is similar to that presented in Sirko & Goodman (2003), except that, while they adopted tabulated opacities (from Iglesias &Rogers 1996 andFerguson 1994), we use analytical opacity fits by Bell & Lin (1994).
The effective temperature is determined from internal heating via viscous dissipation, and scales with : where is the Stefan-Boltzmann constant, Ω = ( / 3 ) 1/2 is the orbital frequency and = (1 − ( min / ) 1/2 ) arises from applying a zero-torque boundary condition at the disc inner radius min . Assuming the disc is in steady state with a constant mass flux, the surface density is given by where = 2 Ω −1 is the kinematic viscosity, assuming that viscosity is limited by disc-scale turbulence, 2 and < 1. The disc is described by a scale height that depends on the sound speed: which also determines the structure and midplane density at each radius by The mid-plane temperature is derived from the diffusion equation, 1 The initial collapse of an optically thin cloud at 0 ∼ 10s K, known as the 'first core' phase, proceeds until the central core reaches H dissociation temperatures (∼ 2000 K). Beyond this stage the protostar continues to grow as a stable object in the 'main accretion' phase (Stahler & Palla 2004). AGN pre-stellar cores are born in the main accretion phase. 2 Turbulence may be driven by magneto-rotational instability, which is expected to occur in a rotating ionized gas and is considered one of the main mechanisms for angular momentum transport in inner AGN discs.
assuming that photons generated in the midplane are transported radiatively to the disc surface according to an effective opacity eff : where eff = 3 8 + 1 2 + 1 4 interpolates between the optically thick and thin regimes, and we have used the Eddington approximation for the photospheric boundary, assuming that the gas is in thermal equilibrium. is the optical depth: computed with the Rosseland mean opacity of the gas at the disc midplane, that is given by a piecewise function following the form = 0 mid (7) where 0 , , and are defined over regimes of and mid following the fits from Bell & Lin (1994). The opacity laws include regimes from ice, dust grains, molecular absorption, Hydrogen recombination, to free-free and bound-free absorption and electron scattering. Opacity scalings and piecewise temperature ranges are summarized in where B is the Boltzmann constant, is the mass of Hydrogen, is the speed of light, the mean molecular weight is taken to be = 0.62 for an ionized gas. For rad we use an effective radiation pressure as in Sirko & Goodman (2003). In the limit of high opacity ( 1), this reduces to the typical scaling with midplane temperature ( rad = 4 4 mid /(3 )). When the gas is optically thin ( < 1), rad ∝ 2 4 mid , and the extra opacity dependence arises as a consequence of the inefficiency of radiation in an optically thin medium.
Once solving the above equations over a range of radii, we can calculate the Toomre stability parameter (Toomre 1964;Goldreich & Lynden-Bell 1965 Beyond the radius where the Toomre becomes less than a critical value 0 ≡ 1.4, the disc is gravitationally unstable. In this regime, we restrict = 0 , effectively holding the disc in a marginally stable, self-regulated state. The solutions for the outer disc carry simple, analytical scalings with chosen parameters, as shown in several previous works (Rafikov 2009(Rafikov , 2015Nayakshin 2006, among others) The assumption of a constant implies a constant midplane temperature and sound speed, and in turn changes the condition for the outer disc density to be = Ω 2 2 0 . (10) Mass conservation means that Eq. 2 still applies, although the source of viscosity mediated by should be no longer interpreted as the scale of turbulent eddies but by the effectiveness of nonaxisymmetries due to GI. The cooling time of the disc is calculated by β crit = 3   Bell & Lin (1994). Given that in regions of low densities multiple conditions can be satisfied, we also set the constraint that > 2 × 10 −4 K for electron scattering.
where we simply take = 5/3 as the adiabatic index. Here we again utilize the opacity at the midplane of the disc, given that our model does not explicitly solve for the vertical disc structure. This estimate for the cooling time comes from the presumption that the energy dissipation due to viscous transport must equal the energy flux through the surface of the disc. Thus, from the equations above, cool ∝ Ω −1 in the inner disc (ensuring steady-state) and drops (in a manner dependent on the opacity) in the outer disc once the assumption of viscous dissipation (via Eq. 1) breaks down with the onset GI.
By assuming the viscosity is dependent on the total pressure, the radiation-pressure dominated region of the inner disc is thermally unstable (Piran 1978;Pringle 1981). One workaround for this is to assume the viscosity depends on the gas pressure only (the so called beta-disc model, Lynden -Bell & Pringle 1974;Huré et al. 2001), which will weaken the viscosity in the inner region and (to maintain the same ) lead to a higher density. For this work we stick with the alpha-disc assumption, given that it does not strongly affect the properties in the outer disc region where gas pressure is dominant.
Full disc solutions are shown for a range of SMBH masses and two values of in Fig. 1. One prominent feature in our disc solutions for lower SMBH masses is an opacity gap at ∼ 10 −3 pc. As we decrease the central SMBH mass, the lower densities and temperatures in the disc allow for the dominant opacity to change from electron scattering and bound-free/free-free transitions to H-recombination, resulting in a sharp opacity drop and a corresponding change in the density profile (further discussed in Section 2.3). Note that the adoption of piecewise opacity powerlaws 3 results in unrealistically sharp transitions in the disc profiles. In a more realistic system, these transitions should be smoother given that opacity regimes may overlap.
We also consider the effect of a temperature background on the disc, which can increase eff and affect the cooling rate. There is observational evidence that surface irradiation contributions a substantial fraction of the disc temperature from reverberation mapping studies of some Seyfert AGN (Fausnaugh et al. 2016;Starkey et al. 2017). If we assume that the temperature background is sourced by the accretion luminosity emitted from the central BH (or the inner disc region), then Eq. 1 obtains an additional component due to external disc heating: Here we assume the disc albedo is zero, meaning that all incident radiation is absorbed, in order to fully demonstrate the effect of irradiation. We find that incorporating this term increases eff by a factor of a few, and has the strongest effect for the BH = 10 6 model. The primary effect is to stabilize the inner disc, or increase the distance to the region where GI and fragmentation set in. Fragmentation occurs where ≤ 0 and cool < crit /Ω. With our disc solutions, we calculate where the conditions for fragmentation are satisfied. Following previous simulations that measure fragmentation conditions, we choose crit = 3. We discuss consequences of the opacity 'gap' for where fragmentation occurs in Section 2.3.

Where does fragmentation occur?
The fragmenting region of the disc is determined by two concurrent conditions, namely that ∼ 0 and that the cooling time is cool Ω < crit , where we choose 0 = 1.4 and crit = 3. It is only when both these conditions are satisfied that a disc can achieve fragmentation. In absence of a short enough cooling time the gas temperature will rise as a result of spiral shocks, raising Q above the fragmentation limit (see Durisen et al. 2007). Beyond the radius where this is satisfied (in a regime Δ frag ), the disc will fragment into bound clumps which accrete surrounding material as they collapse, forming protostars. Note that the disc does not necessarily need to be more massive than the central SMBH in order to become self-gravitating: if it is thin and dense, the criteria are satisfied before disc ∼ 10 −2 BH (see Fig. 2). While theory predicts that discs around more massive BHs are more for BH = 10 6 M and = 0.1 for accretion rates at different fractions of the Eddington rate ( / Edd = Edd ). Star symbols delineate the transition radius t , or the innermost region where fragmentation occurs. In the inner disc, opacity is dominated by (from left to right) electron scattering, then bound-free transitions until a gap occurs due to H − absorption, and then molecular lines and dust scattering. Fragmentation occurs within this gap or beyond it, depending on the accretion rate which sets the disc temperature. prone to GI (Goodman 2003;Lodato 2007), we find that discs around less massive SMBHs may fragment at smaller radii, provided the accretion rate is sufficiently high ( Edd > 0.05). This is not intuitivedespite the lower gas densities, the cooler temperature allows for a different opacity source to dominate and consequently a faster cooling time. Note that if you consider the fragmentation radius in terms of the respective Schwarzschild radii of the central SMBH, for which S ∝ BH , then more massive SMBHs fragment at distances of fewer Schwarzschild radii. However, the lower densities in discs around less massive SMBHs mean that the fragments (and subsequent AGN stars) are born with comparatively smaller masses, as we will show in Fig. 5. The inclusion of surface irradiation increases the distance to the transition radius by a small factor. For the parameter space explored in this work, it does not prevent fragmentation. . Transition radius t from inner to outer disc (or minimum radius where = 0 and Ω cool < crit ) for discs with = 0.1 around a range of SMBH masses. Different lines/symbols correspond to solutions for the corresponding constant accretion rate = Edd Edd . The hollow circles correspond to the Edd = 0.1 disc with surface irradiation, which pushes fragmentation out to farther radii. Due to opacity transitions at lower temperatures, discs around 10 6 M SMBHs can fragment at smaller distances depending on the accretion rate.

Centi-parsec fragmentation near opacity gaps
The transition to the fragmenting region of the disc and whether or not this occurs within (or due to) the opacity gap depends on the disc accretion rate. Fig. 3 shows solutions for surface density and opacity profiles for discs with = 0.1 around a BH of mass BH = 10 6 M , for three different Eddington ratios. Discs at higher steady-state accretion rates are denser, leading to GI and fragmentation at closer radii ( t ∼ 10 −3 pc, denoted by the star symbols). The solution with the lowest accretion rate ( Edd = 0.05) best demonstrates the full opacity gap, which is comparable to solutions found in Thompson et al. (2005). A lower accretion rate leads to fragmentation at farther radii ( t 10 −2 pc) beyond the opacity gap, in a region dominated by dust scattering and sublimation. In the bottom panel of Fig. 3 we show the opacity profiles for each disc model, where the star symbols indicate beyond which radius the disc is gravitationally unstable. In the inner disc, the opacity is dominated by electron scattering and bound-free+free-free transitions. At intermediate radii, temperatures become low enough for H − recombination to occur, which results in a sharp opacity drop, beyond which molecular line cooling begins to dominate. At cooler temperatures, the opacity becomes dominated by dust sublimation and dust scattering. In these regions the disc is still optically thick, as shown in Fig. 1, and thus the transition in opacity strongly affects the cooling rate. In summary, the dominant opacity law within the fragmentation region may vary from H-absorption to molecular interactions to dust grain sublimation, depending on the combination of density and temperature at t . 4 In Fig. 4 we show the minimum radius of fragmentation, or the transition radius t , versus the central SMBH mass for different accretion rates. Differences in opacity laws determine the cooling rates and lead to the scatter at lower SMBH masses. For all solutions shown here, fragmentation begins to occur within 10 −3 10 −1 pc, but, as we show in the following Section, the variance in disc properties in this region lead to different initial clump properties.

INITIAL MASS DISTRIBUTION OF AGN STARS
In this section we describe how the properties of the disc in the gravitationally unstable region lead to predictions for the size of clump formation, accretion rates onto pre-stellar cores, and an initial mass distribution of AGN stars.

Fragmentation into bound gas clumps
In regions where the disc fragments, gravitationally bound clumps 5 arise with a range of sizes centered around a characteristic scale. In the linear theory, this characteristic scale mu is derived from the dispersion relation for local axisymmetric density perturbations, and corresponds to the most unstable wavelength The maximum scale at which fragmentation can occur is given by twice the most unstable wavelength. Pre-stellar clumps are thus typically born with characteristic initial radii c = 1 2 mu , and the corresponding enclosed mass 6 is T = Σ where we assume a circular volume element with surface density Σ, and we have plugged in a normalisation corresponding to the outer disc properties for a 10 6 M SMBH. As discussed in Section 2.1, the conventional theory over-predicts the size and mass of clumps when compared to simulations of magnetized or hydrodynamic selfgravitating discs. We include it here as a 'maximum estimate' in order to compare to previous work. An improved estimate that takes into account the fact that, in an differentially rotating disc, the instability has a global rather than local character and occurs in non-axisymmetric conditions, namely along spiral arms (Boley et al. 2010) finds that the perturbations are dependent on the density enhancement in the spiral arm, which corresponds to a scale given by T = T /M 2 assuming the flow in the arm behaves as an isothermal shock of Mach number M. With this estimate, and further assuming finite thickness across the arm, ( J /2) 3 where J = ( 2 /( )) 1/2 , this further suggests that the initial clumps formed by GI are already too large to maintain stability, and hence should undergo collapse. However, the Jeans criteria is oversimplified for star formation in the ISM, as it underpredicts the characteristic mass when compared to the observed mass function (e.g.
where g is a shape factor of order unity. We consider frag a fiducial estimate for the initial mass of prestellar fragments, which notably predicts masses that are an order of magnitude less massive than . We show the range of clump masses as a function of the fragmentation radius for each disc model in Fig. 5. Note that each line represents a value of the characteristic clump mass at each radius, rather than a strict prediction of the mass distribution. A few to several clumps may form in each radial bin with a distribution centered on the characteristic mass. The initial mass distribution of the population will fall within the given mass range, but with a slope that depends on the gas density and SF efficiency.

From clumps to protostars
Once bound, pre-stellar clumps may continue to collapse depending on the gas properties. If the central temperature of a fragment is below the Hydrogen dissociation temperature (∼ 2000 K), collapse of a fragment is initially quasi-static and slow (Decampli & Cameron 1979). Once molecular Hydrogen dissociates, the collapse of fragments occurs dynamically (Larson 1969). The former contraction phase is a critical step in conventional star formation (as well as for planet formation in protoplanetary discs, Helled et al. 2014), during which the first core forms. In this phase, the surrounding gas is vulnerable to internal feedback once the cloud is optically thick, which slows down the infall of surrounding gas. In the AGN case, however, star formation occurs in an environment with initially high temperatures (see the . Components of the pre-stellar clump in a non-rotating, isolated case (top) and a disc-embedded, rotationally-supported case (bottom). The clump is a gravitationally bound, collapsing cloud. Neglecting rotational support, the protostar embryo consists of a central hydrostatic core and a surrounding free-falling envelope. For the rotationally-dominated case, the mass of the clump is partly in a central core (the nascent protostar) and partly in a surrounding disc. Additional gas can be supplied to the circumclump disc via gas inflow into the Hill sphere. Accounting for rotational support results in initially less massive stars. midplane temperature in Fig. 1 7 ), and critically, the bound clouds are born in a regime where thermal effects from their contraction are unimportant. Instead protostars are effectively born in the main accretion phase, which is analogous to the dynamical collapse phase in planet formation. However, we stress that this outcome is where the analogy between fragmentation in protoplanetary discs and in AGN discs breaks down -in the planet case, protoplanets are vulnerable to disruption during their initially slow contraction (Müller et al. 2018), while in the AGN case, the initially high temperatures and 7 In principle, AGN stars may form in gas that is below the H dissociation temperature, which can occur in discs around less massive SMBHs (in fact, the 10 6 M , = 0.1 disc model in Fig. 1 shows that mid ( frag ) < 2000 K).
In such cases the collapse may be initially delayed but should quickly presume dynamical collapse.
densities ensure that clumps collapse dynamically to form stars and are unimpeded by internal feedback or subsequent disruption. During (spherical) dynamical collapse, a gravitationally bound clump will collapse on the free-fall timescale: ff = 3 32 1/2 ≈ 10 0.01pc where is the average density of the cloud, and for the right hand side we have plugged in the disc density scaling from Eq. 10, which provides an upper limit of the free-fall time (since the protostellar clump is ab initio denser than the disc). The estimate is derived for the time for a cloud to collapse to infinite central density, but inevitably before this point the inner temperatures will become sufficiently high to ignite nuclear fusion. The rapid free-fall timescale, a result of the high densities in the disc model, suggests that clumps in the fragmentation region collapse swiftly into AGN stars within 10 0 ff 10 3 yr, depending on the initial density of the gas. An additional consequence of such rapid collapse is that fragments become relatively compact, reducing their interaction cross section soon after formation. This suggests that it is unlikely for fragments to merge with each other during this stage of formation, unless a relatively large number of fragments form.

Accretion during protostar collapse
During the main accretion phase, protostars continue to grow via accretion from their surrounding envelopes. To provide context on the scale of rapid growth, we describe two scenarios for this process: the first a simplified, inside-out spherical collapse (considered a maximum estimate) and the second is an axisymmetric collapse in which the growth of the prestellar core (PSC) is governed by the rate of angular momentum dissipation through a surrounding accretion disc. The latter case further supports the assumptions we adopt in this work, in which the initial clump masses are lower than previously anticipated. In Fig. 6 we provide a schematic of the pre-stellar clump assuming purely spherical infall (top panel) and including a rotationally supported accretion disc (bottom panel). The latter is a natural expectation from fragmentation in an accretion disc. In a shearing medium, PSCs are inevitably born with angular momentum. Higher angular momentum material settles into a protostellar disc, which can delay collapse and mediate the growth of the core (Galvagni et al. 2012).
In the spherically symmetric case, a PSC accretes from an envelope at a rate governed by the local sound speed, also known as 'inside-out collapse' (Stahler & Palla 2004 where is the sound speed of the gas in the co-orbital region and on the right side we have plugged in the parameter scaling implied by = 0 , assuming the central BH accretes at a constant fraction of the Eddington rate. This expression assumes the envelope is isothermal and non-rotating, but the infall rate has a weak dependence on rotation (Stahler & Palla 2004). Given the high temperatures of the AGN disc, inside-out collapse predicts extremely rapid accretion onto protostars ( 0.1−10 M yr −1 , depending on BH ). However, this accretion phase is short-lived and only occurs within the sphere of influence of the core, until the nearby gas is quickly depleted.
The growth of prestellar cores that have significant rotation proceeds differently, yet may lead to similarly high accretion rates. In this case accretion occurs via a circumstellar disc, which regulates the flow onto the protostar at a rate set by viscous dissipation (the details of which we neglect here) as well as the gas inflow from the surrounding environment. If we consider the shear flow around the protostar, accretion into the Hill radius occurs at a rate (following Boley et al. 2010) where we have assumed that the surface density is constant within the protostar region and normalised to values of interest to our system. Not necessarily all of this mass is accreted by the star, given that gas can flow in and out of the Hill radius. Furthermore, accretion via a disc must be regulated by viscous dissipation, which is often inefficient (Shabram & Boley 2013), so Eq. 18 should be interpreted as an upper limit. For both the spherical and axisymmetric cases, accretion rates onto collapsing cores are determined by gas properties which, as a consequence of our disc model and the scaling of frag , are constant for all nascent stars.
For perspective, note that these accretion prescriptions can exceed the Eddington rate for an object of mass * , Note that we use the electron scattering opacity here, but in cooler regions of the disc where dust can condense, the lower opacity (as shown in Fig. 3) leads to a correspondingly higher possible Eddington rate.
The accretion rate during this phase is high but short-lived, given the rapid free-fall timescale (indeed a strong limit also arises from the local gas supply which quickly runs out). Protostars only increase their masses by at most a factor of ∼ 1/10th during collapse (assuming they accrete over a timescale ff ), resulting in characteristic stellar masses within the range determined by the initial clump mass (Eq.15), as shown in Fig.4. 8 In principle, gas inflow from the global accretion disc could continue to feed the circumstellar disc after star formation. At this point we assume that once the star is formed, feedback (either from accretion or stellar winds) would suppress further growth, but we discuss a possible outcome of subsequent stellar accretion in Section 3.4.

Opacity-limited mass
A minimum clump mass limit can be calculated by determining the point at which the clump becomes opaque to its own radiation. This minimum mass is known as the 'opacity limit' mass (Low & Lynden-Bell 1976;Rees 1976;Silk 1977). Assuming spherical collapse, the opacity-limit mass is of order min = 0.025(1/ ) 16/7 ( / ES ) 1/7 M , where is the mean molecular weight and the final opacity at which the fragment becomes opaque. The limit is weakly dependent on the opacity and remains 8 In the case of SMSs, the accretion rate in the main accretion phase is important for determining the final stellar mass (Hosokawa et al. 2013) and the stability of the growing core (Haemmerlé et al. 2019). By comparing to core accretion rates of models by Haemmerlé et al. (2019), one can see that even these rapid growth rates remain in the regime of 'hydrostatic growth,' i.e. the structure of the PSC can readjust as it gains mass and avoid runaway growth that leads to instabilities. below min 10 −1 M for 1cm 2 g −1 , which is consistent with the range of fragment masses predicted by Eqs. 14 and 15.
We note that this limit may change if one considers that fragmentation occurs in the vicinity of an accreting SMBH which can impart a temperature background (if the fragmenting material is optically thin). In the presence of a radiation field, the opacity-limited mass becomes more sensitive to the opacity under the assumption that it cannot cool below the background temperature (Low & Lynden-Bell 1976). We defer this calculation to future work, noting that the precise opacity-limit mass in the AGN regime is a complex radiative transfer problem that depends on the composition of the surrounding material, the adiabatic evolution of the fragment, and the preferential geometry of the emitted feedback, from both the AGN and the collapsing fragments.

From disc properties to the stellar mass distribution
While the mass distribution shown in Fig. 5 demonstrates the range of possible stellar masses during a fragmentation episode, it does not quantify the total number of stars formed. For this we must assume a SF efficiency , which for simplicity we assume is constant with radius. We adopt = 1% (or = 0.01), which is motivated by efficiencies calculated in models of starburst discs (Thompson et al. 2005) 9 and supported by observations of starburst galaxies (Kennicutt 1998) and giant molecular clouds (GMCs) Murray (2011). (We caution, however, that the star formation considered in our model occurs on smaller spatial scales, and need not necessarily lead to similar efficiencies as observed in GMCs.) Similarly, the efficiency of planet formation via GI in protoplanetary disc simulations suggests 10% of the gas goes into bound clumps, though the fraction is lower if one accounts for subsequent disruption or additional nonlinear effects (Galvagni et al. 2012;Müller et al. 2018;Schib et al. 2021). To remain conservative (and to ensure that star formation does not shut off the accretion flow), we adopt the fiducial value of = 0.01. To account for the possibility that star formation is more efficient, e.g. in the case that feedback from star formation is not sufficient to restabilize the disc, we also consider the case of = 0.3, motivated by the maximum star formation rates observed in GMCs. This choice has implications for the final EMRI rate, which we discuss in Section 5.
Stars form over a timescale proportional to the orbital time: ff ∼ Ω −1 , so the efficiency provides a star formation rate * = 2 where Σ * = ΣΩ is the star formation rate per unit area, integrated over the fragmenting region from t to the outer disc edge out . The efficiency can be interpreted as a reduction in the timescale of star formation, or, equivalently, tells us what fraction of the gas turns into stars. If no additional gas is supplied to the disc, fragmentation will result in a population of * stars * = out t where / is the combined mass of stars formed within a radial bin , which is dependent on the gas density: The differential mass distribution / also depends on radius. In each radial bin, stars form at masses sampled from a distribution: where ( , ) is the radially-dependent mass spectrum for which we choose a normal distribution where ( ) is the standard deviation, taken to be = frag ( )/3 (to ensure that the masses do not exceed twice the characteristic mass), and = frag ( ) is the distribution mean. This results in a combined mass of stars * ,tot = out t = 2 and the radial number distribution is constrained by the mass budget as follows: Stars in each radial bin / follow a mass distribution / centered on the characteristic mass frag ( ). In practice, the distribution is obtained by splitting the radial profiles shown in Fig. 4 into 50 discrete, logarithmically-spaced bins, and sampling ∼ / frag stars from a distribution centered on the characteristic mass frag ( ) within each bin. The final distribution does not depend strongly on the number of bins, as long as enough gas is contained within each to form a finite number of stars. We assume stars form only as single systems, neglecting the binary fraction.
This results in the initial mass function shown in Fig. 7, for which we plot models with = 0.01. The AGN stellar mass distribution is top-heavy compared to the conventional Salpeter initial mass function (IMF): Fitting a line to the data shown in Fig. 7 provides a slope / ∝ Γ , with Γ ∼ −0.7. This is considerably flatter than the Salpeter IMF, where the slope is Γ SP = −2.35. Our finding is in agreement with the stellar mass distributions measured in simulations by Mapelli et al. (2012) of a fragmenting disc for a Sgr A* type SMBH. For different gas masses and thermodynamic assumptions, they measure an IMF slope ranging from Γ ∼ −0.5 to −1.6. In another related set of simulations by Nayakshin et al. (2007), fragmentation of a disc around Sgr A* also produces a top-heavy IMF that is dependent on the stars subsequent accretion of the surrounding gas. Such a top-heavy IMF is also supported by observations of ubiquitous metal enrichment in AGN, which can be produced by supernovae of massive stars (Wang et al. 2011;Toyouchi et al. 2022;Xu et al. 2018;Lai et al. 2022;Dittmann et al. 2022).
Our choice of a constant SF efficiency (e.g. Eq 22) requires that the number of stars formed decreases at outer radii to ensure that the gas disc is not depleted. The consequence of this is that the star forming region of the disc is finite, and only a small number of stars form at outer radii, where the heaviest end of the distribution lies. Note that in the case of a higher SF efficiency, the amplitude of the distribution increases according to Eq. 22, while the the range and slope of the distribution remains the same. In principle the disk could extend to larger radii, in which case a higher SF efficiency would potentially increasing the heavy end of the distribution and the resulting EMRI rate. To keep the model conservative, we choose to set the outer disk edge to be out = 0.1 pc, which is also in agreement with constraints on disk sizes from observed spectral energy distributions (see Sirko & Goodman 2003). In the case that the efficiency varies with radius, the slope of the IMF will change (as will the steady-state profile of the disc). Furthermore, we neglect here the possibility that fragments may merge, an increasing possibility for cases of higher , which would skew the distribution to higher masses.

Bondi accretion onto embedded stars
Here we demonstrate that despite the lower masses of in-situ stars, accretion from their environment at limited rates can rapidly skew the mass distribution to higher values. For stars moving on prograde, circular orbits (such that the relative velocity between the orbiter and the gas is negligible), the Bondi-Hoyle accretion rate is given by (Hoyle & Lyttleton 1939;Bondi & Hoyle 1944;Edgar 2004 where B = 2 * / 2 is the accretion radius, which is defined by the region of gravitational influence. This estimate suggests rapid accretion onto embedded stars, and has been presented as a rapid growth mechanism for gas-embedded stellar populations (Davies & Lin 2020). In fact, the disc properties derived above suggest a Bondi accretion rate that exceeds the accretion rate onto the central SMBH, which would break more than the steady state disc assumption. However, the Bondi radius is in most cases larger than the disc thickness, and thus the assumption of rapid, spherically symmetric infall is not appropriate, especially if stars grow to larger masses. While accretion inside the Bondi radius may occur in a super-Eddington configuration-certain accretion/feedback geometries allow for super-Eddington growth, e.g. a slim disc (Abramowicz et al. 1988) or accretion with non-equatorial feedback (Jiang et al. 2014) or reduced radiative efficiency (Jiang et al. 2019)-the rate must at least be limited by the amount of gas that can flow into (and stay within) the star's sphere of influence.
To demonstrate an outcome of stellar accretion, we apply a modified Bondi rate to the initial population for which the accretion radius is limited to either the star's Hill radius or the disc thickness, such that acc = 1 2 minimum( Hill , ℎ( * )) where Hill = * ( * /3 BH ) 1/3 , ℎ( * ) is the disc scale height at the stars position * , and the factor of 1/2 ensures that the accretion radius is well within the disc 10 (and consequently so that the rate remains at least an order of magnitude below the disc accretion rate). The accretion rate is given by where we introduce an efficiency factor of = 0.01 to account for empirical evidence that gas flowing into the sphere of influence does not guarantee its accretion onto the star. In fact, simulations show that gas flows past an embedded migrator from the outer to the inner disc, (see e.g. Duffell et al. 2014or Derdzinski et al. 2019). Fig. 8 shows the evolution of the stellar mass distribution after ∼ 10 7 yr, a timescale set by the AGN disc lifetime discussed in Section 4. Here the accretion rates are determined by the disc properties at the initial position and updated as the stellar mass increases. This assumes the lowest disc densities, leading to a conservative estimate of the Bondi rate, but also the highest disc thickness, which allows for a larger accretion radius. This accretion increases stellar masses by more than two orders of magnitude throughout the disc lifetime. The result is more pronounced for discs with = 0.01, for which the higher disc density leads to more rapid growth (in these cases we limit the stellar masses to 10 4 ). The shift in peak of the stellar mass distributions for all models is shown in Table 2.
We consider this accretion rate an overestimate, given that during growth, stars will likely reach a critical mass beyond which accretion becomes inefficient due to radiative feedback or momentum-driven winds. The precise limit depends on details of stellar winds, feedback, and the accretion geometry, but may lie between 30 to 60M (Ekström et al. 2012;Cantiello et al. 2021, although see Dittmann et al. 2022). Due to the uncertainties regarding stellar accretion, in the following section we use the initial mass distribution (shown in Fig. 7) as our fiducial population, which is based on the fragment mass limits from Eq. 15. As we show, including stellar accretion will alter the final numbers when determining observational outcomes, but the basic story remains the same: stars of different masses will have different fates due to the interplay of migration, stellar evolution, and tidal disruption.

SUBSEQUENT EVOLUTION AND MIGRATION OF AGN STARS
We now turn to the subsequent migration of AGN stars. PSCs collapse on a timescale faster than migration (see the dynamical time in Fig. 9), which allows us to calculate the subsequent evolution of AGN stars separately from their formation. We break down the range of outcomes for AGN stars into three categories: (i) Tidal 'disruption': Stars migrate inwards until they reach the tidal radius where they are destroyed.
(ii) In-situ EMRIs: Stars migrate inwards until they are accreted by the SMBH, either by avoiding disruption (depending on the relationship between the SMBH mass and stellar structure) or collapsing to compact remnants.
(iii) Leftover stellar population: Stars migrate inward until the disc disperses, surviving the AGN phase.
Each of these outcomes is dependent on the mass of the star, the mass of the SMBH, and the disc properties. For the present work we neglect details that would occur as sub-categories (e.g. whether an EMRI occurs as a star or a specific stellar remnant, it is categorized as outcome (ii)). Our fiducial model assumes that no further accretion occurs after star formation, meaning that we adopt initial mass distributions in Fig. 7 and keep the stellar masses constant. We discuss implications of including stellar accretion (as shown in Sec. 3.4) in Sec. 6.
Stars embedded in the disc will excite nonaxisymmetric perturbations that react on the star's orbit, resulting in an exchange of angular momentum that will affect the orbital separation (Goldreich & Tremaine 1980). In the limit of small mass ratio (here neglecting dependencies on gas thermodynamics), this migration regime is referred to as 'Type I,' and the strength of the (typically inward) torque on an embedded migrator of mass * at a separation in the isothermal limit can be estimated with the following expression (Tanaka et al. 2002): where I is a factor that depends on disc gradients. Rather than solving for I from our disc solutions which can lead to sign changes in the torque, here we fix I = 2 to represent a typical inward migration rate. In Appendix A, we discuss the implications of disc gradients for migration reversal. This corresponds to a change in radial separation I = 2 I / / * , where = ( BH / ) 1/2 is the Keplerian orbital velocity.
In the classical paradigm, migration slows down if the migrator carves a gap, reducing the density in the coorbital region, which weakens the migration torque. In our case, gap-opening for AGN stars is unlikely for stars with small mass ratios compared to the central BH (e.g. for 10 −5 ), as suggested by conventional gap-opening criteria (Crida et al. 2006). Furthermore, a number of numerical studies conducted for self-gravitating discs on various scales have shown that migration in such systems can be even faster than implied by Eq. 30 due to the contribution of torques from asymmetric flows in the vicinity of the perturber (Malik et al. 2015;Mayer 2013), which either inhibits gap formation completely or at least suppress the formation of a deep gap that can slow down migration sensibly (Szulágyi et al. 2017). We note that the migration of a population of stars may not be well-represented by expressions that assume a single migrator. Simulations tend to show that while a fraction of fragments tend to migrate inward at a rate comparable to the Type I rate, other fragments may experience migration that slows down or reverses direction (Helled et al. 2014). In the inner disc, gap-opening or partial gap-opening may occur, especially if stars are able to accrete and grow. For simplicity we consider only inward migration but with two limiting estimates (Type I and a slower, gap-corrected expression).
In the case that a gap begins to open, numerical simulations find in many cases that the torque on a perturber scales with the reduced surface density in the gap, and can be estimated by (Fung et al. 2014;Kanagawa et al. 2018) where = 2 ( /ℎ) 5 −1 . This estimate takes into account the possibility of partial gap-opening, allowing for a continuous expression for the torque as a function of and disc parameters. When calculating migration rates for the stellar population, we take into account the reduction in gas density due to star formation. This only results in a minor effect for models with high SF efficiency ( = 0.3), for which the migration times increase by a small factor 2. We caution that the above expression assumes a non-self gravitating, isothermal accretion disc, and is not tested for all disc properties or thermodynamics. Furthermore, it is known that gas dynamics in the gap-opening regime become highly nonlinear, and the resulting torque can have nontrivial dependence on disc parameters (Duffell et al. 2014), especially when resolving gas close to the perturber (Derdzinski et al. 2019) or when gaps are particularly steep (Chen et al. 2020). However, due to the majority of stellar masses being below the gap-opening threshold, the correction has minor implications for the outcome of AGN stars (discussed in Section 5).
Additionally, orbital angular momentum loss occurs due to GW radiation from interaction with the central SMBH, which, assuming a circular inspiral, results in a separation evolution (Peters 1964) where is the gravitational constant, and is the speed of light.
The timescale for a star to migrate from its initial separation i to a final separation f is mig = f i 1 I + GW (33) assuming that gas-induced migration and GW migration are uncoupled. More massive stars migrate faster and reach the inner disc first, despite being born at farther separations. In theory migration should slow down once the star becomes a less massive compact remnant. We neglect this detail here as it will only affect when, not if, the stars produce EMRIs.
Stars will migrate till they are accreted by the central SMBH or disrupted by tidal forces. Thus the outcome of stars relies on the interplay between the migration timescale and the stellar lifetime. Tidal disruption 11 occurs if stars reach the tidal radius, or the region where tidal forces from the central SMBH overcome the gravitational binding energy of an incoming star, defined by where * is the radius of a star with mass * (Hills 1975). The expression is approximate due to uncertainties in stellar structure, for which we estimate stellar radii by * = * M 1/2 R (35) 11 Tidal disruption typically requires that stars reach the tidal radius on an eccentric orbit (which is also a consequence of dynamical formation scenarios in the absence of gas), and in this case the star is violently disrupted. In the case of disc-embedded stars with near circular orbits, it is more likely that stars are tidally dissolved, or experience Roche lobe overflow (Rossi et al. 2021;Dai & Blandford 2013;, and that torques from the disc may complicate the picture. We can conveniently still refer to this event as a tidal dissolution event (TDE).
which underestimates the size of sub-solar mass stars 12 as well as massive stars, at least in isolated environments (see e.g. Hosokawa et al. 2013). We neglect partial disruptions, which may occur for massive stars with extended envelopes (Rossi et al. 2021;MacLeod et al. 2012) and tidal heating, which may inflate stellar radii (Alexander & Morris 2003;Li & Loeb 2013). Tidal disruption is the primary obstacle for survival of stars around less massive SMBHs ( BH 10 7 M ), for which the tidal radius lies outside of the SMBH event horizon. Disruption can be avoided if stars collapse to become CRs before migrating to this regime. The stellar lifetime also depends on the stellar mass: where follows from the mass luminosity relation ( * ∝ 1− * ). For stars below * /M < 15, we adopt the main sequence values = 2.5 and = 10. For stars more massive than * /M ≥ 15, we adopt scalings from rotating massive stellar models at solar metallicity (Ekström et al. 2012, but see also Yusof et al. 2013) for which = −0.78 and = 8.00.
The final timescale for comparison is the AGN disc lifetime, which is given by the viscous time at out , assuming that no additional gas is supplied to the disc during the accretion episode: where we have plugged in the turbulence assumption for kinematic viscosity.
The migration times shown in Fig. 9 assume a constant stellar mass, without taking into account that beyond a stellar lifetime, stars 'become' a compact remnant of smaller mass. This does not change the migration time considerably, since most of the migration occurs during the stellar lifetime. The main implication of the remnant mass is whether the EMRI occurs during the AGN disc lifetime or after disc dispersal. More massive remnants will inspiral more quickly due to stronger GW emission, and thus more likely occur while the disc is still present.
In Fig. 9, we show the timescales for Type I migration (dark purple lines) from the starting separation to either the tidal radii or to 10 S , where S = 2 BH / 2 is the Schwarzschild radius, depending on which occurs first, as a function of the stellar mass. For comparison we show the dynamical time at the formation radius (grey line), the stellar lifetime (green line), and the AGN disc lifetime (blue dashed line). Intersections in the lines determine the minimum mass below which stars do not survive to leave CRs or form EMRIs. Migration is also shown with the gap correction (faded purple lines, using Eq 31), which slows down migration for stars with mass ratios greater than 10 −6 , suggesting that for each disc model more massive stars survive disruption.
Gas-driven migration on its own becomes inefficient at inner radii (at 10 −5 pc), where the disc is thick and the gas density is low. However, in this regime GWs become important and soon dominate the inward migration, leading to comparatively fast migration timescales, around mig ∼ 10 6 − 10 7 yr depending on the disc model and stellar/remnant mass. In most cases, the migration time to disruption or accretion is less than the stellar lifetime. In the case of disruption, this suggests that stars do not survive to leave CRs or produce in-situ EMRIs.
12 Sub-stellar objects or brown dwarfs may survive disruption, contributing to the 'eXtremely small mass ratio inspiral event rate (see Fig. 1 in Amaro-Seoane 2019

Other migration regimes
Simulations of Type I migration in the non-isothermal and isothermal limits find more detailed expressions for the torque that depend on the density, temperature, and entropy gradients of the disc (Kley & Crida 2008;Tanaka et al. 2002;Paardekooper et al. 2010). With our disc solutions we can calculate the torque profiles and determine for which regions the isothermal or adiabatic torque regimes occur. These calculations and resulting torque profiles are shown in Section A of the Appendix. We find that the magnitude of torques is within an order of magnitude of the Type I estimate from Eq. 30, whether it is in the adiabatic or isothermal regime. More interestingly, our disc solutions suggest the presence of regions where the torque changes direction from inward (reducing angular momentum of the migrator) to outward (increasing its angular momentum), typically at distances ∼ 10 − 10 2 S from the SMBH. The intersection of these regions is referred to as a 'migration trap' in the context of both planets (Lyra et al. 2010) and BHs (Bellovary et al. 2016).
If star and/or compact object migration is halted in these regions, this could increase the interaction rate between migrators (Secunda et al. 2019). However, it is likely that migration continues due to other processes and that such traps are not sustained (e.g. GW emission). In particular, the highly dynamic regime of self-gravitating discs is unlikely to preserve any feature requiring a delicate torque balance. Motivated by the latter, we neglect migration traps in our calculation of the stellar orbital evolution and resulting EMRI rate. Instead we adopt the assumption that migration is fast, inward, and rapid, following the scaling in Eq. 30, which is a choice more in line, qualitatively, with how migration is currently understood in self-gravitating discs. Type I migration leads to a conservative in-situ EMRI rate, as the inclusion of the gap correction or traps will delay the migration of stars, providing more time for them to collapse to compact remnants. It is also supported by recurring evidence of fast migration in discs with GI (e.g. Malik et al. 2015, but see the code comparison paper by Fletcher et al. 2019).

OUTCOME OF THE AGN STAR POPULATION
By comparing the migration timescales to stellar lifetimes and AGN disc lifetimes, we can categorize outcomes of AGN stars. We focus primarily on in-situ EMRIs before briefly touching on TDEs and leftover stellar populations.

In-situ gas-embedded EMRIs
Under the assumption that no additional gas is fed to the disc, we calculate a minimum EMRI rate per AGN during an isolated accretion episode. Such a scenario is plausible if accretion onto the SMBH produces feedback that inhibits additional gas from falling into the SMBH Bondi radius. The rate is an underestimate if SMBHs experience multiple episodes of accretion over longer timescales.
The in-situ EMRI rate depends on the efficiency of star formation (how much of the gas turns into stars within Δ frag ), the stellar mass distribution (determined by disc properties), and subsequent migration rates. From one star formation event, we can quantify the total number of stars formed as * (Eq. 21). Whether stars become EMRIs is determined by the timescale comparisons in Fig. 9. In order to be accreted by the SMBH within the disc lifetime, stars must satisfy three constraints that depend on their masses and initial separations: (i) Stars migrate toward the SMBH with the AGN disc lifetime Figure 9. Migration timescale from intial separation to 10 (solid purple lines), or to tidal (dashed purple lines), stellar evolution timescale (green line), AGN lifetime (blue dashed line), and orbital time (grey line) for three different SMBH masses and two disc viscosities. Migration times are shown for Type I (dark purple line) and with a gap-opening correction (faded purple line) which significantly slows down migration for higher stellar masses. Stars in the left panels are disrupted before they leave compact remnants, while higher mass stars in the middle and rightmost panels survive disruption to be accreted, producing in-situ EMRIs.
( mig < AGN ), which means that their mass must be above a critical mass i .
(ii) Stars avoid disruption because the tidal radius (Eq. 34) is within the regime where they become GW sources ( tidal < 10 S ( BH )). This translates to their mass being below a limit ii (e.g. as shown in Fig. 9 for BH = 10 7 M ).
(iii) Stars avoid disruption because they collapse to compact remnants before reaching the tidal radius ( * < mig ), which requires that their masses are greater than a limit iii so that their lifetime is more rapid than the migration time.
In other words, EMRIs originate from stars where the migration timescale (solid purple lines in Fig. 9) is less than the AGN disc lifetime, excluding regions where stars are disrupted (dashed purple lines). The total number of EMRIs per AGN is given by summing up the stars from the mass distributions (e.g. Fig. 7) within the mass range that satisfies these conditions, such that emri = * [( * > i ) and (( * < ii ) or ( * > iii ))]. Over an accretion disc lifetime AGN , the in-situ EMRI rate per accretion episode is We show this EMRI rate for all disc models in Fig. 10 and in Table 2. Not all models produce in-situ EMRIs. Specifically, the survival of stars in the 10 6 M BH disc model is sensitive to assumptions of stellar accretion and migration. Stars formed in this system are at lower masses, and efficient (Type I) migration leads to their disruption within ∼ 10 6 yr. If stars can accrete to larger masses, however, or if migration slows down due to partial gap-opening, this allows for a subset of the population to produce compact remnants. In these cases, the EMRI rate is emri ∼ 10 −7 −10 −5 yr −1 . For disc models with larger SMBH masses, in-situ EMRIs occur at rates emri ∼ 10 −5 −10 −4 yr −1 depending on the disc properties, the torque prescription, and the SF efficiency. The inclusion of stellar accretion boosts the EMRI rate, in some cases by an order of magnitude. The effect of the gap-correction is less clear: in principle it allows for more massive stars to survive disruption and produce EMRIs (this is seen in the bottom middle panel of Fig. 9). In practice, most stars are below the mass where the gap-correction becomes relevant, as shown in Sec. 3, unless they accrete substantially. The higher viscosity disc (with = 0.1) for BH = 10 7 produces more EMRIs, despite the stars having slower migration times. This is due to the lower density discs producing larger numbers of less massive stars that can avoid disruption. In a similar fashion, models with high SF efficiency lead to a higher EMRI rate, simply by increasing the number of massive stars that survive to seed EMRIs.
We note that by limiting the disc size to 0.1 pc, we also limit the number of massive stars that form, and in turn limit the EMRI rate. We expect these rates will increase for larger discs that experience higher SF efficiencies. On the other hand, these rates may decrease if one considers mutual gravitational scattering between stars/BHs, which may preserve some stars from disruption but also lose potential EMRIs (Forgan et al. 2018). Interactions between BHs which can merge with each other via gas-assisted binary formation (Tagawa et al. 2020a) could also decrease the rate, while producing higher-mass EMRIs (or IMRIs).
The rates we find are in agreement with the rate of in-situ EMRIs predicted by the models of Dittmann & Miller (2020). For further comparison, , compute EMRI rates due to interactions of nuclear cluster stellar/BH population with a massive accretion disc. For the 'conventional' EMRI scenario in dry nuclei (via two-body interactions), they find EMRI rates around ∼ 10 −6 [yr −1 ] per galactic nuclei. In the presence of a gas disc, the EMRI rate typically increases to ∼ 10 −6 − 10 −5 yr −1 per AGN, with an order of magnitude variance depending on the mass of the SMBH and the model parameters. These numbers are in line with our calculations for models that produce EMRIs, suggesting that the in-situ channel can generate EMRIs at a comparable rate to the orbit-capture channel. Tagawa et al. (2020a) consider the evolution of a population of stellar cluster stars and BHs that are captured into binaries in an AGN disc. This work is focused on the evolution of stellar mass binaries via gas-assisted mergers, but also predicts a consequent EMRI rate 0.1 − 0.6 Gpc −3 yr −1 . For comparison, we can express the EMRI rate in terms of the central SMBH accretion rate: m = emri / BH ( AGN ) M −1 . This can be used to derive a volumetric in-situ EMRI rate˜= m SMBHs , where SMBHs is an estimate of the total accretion onto all SMBHs in the local universe. Calculations by Marconi et al. (2004) suggest SMBHs ∼ 5 × 10 −6 − 5 × 10 −5 M yr −1 Gpc −3 for ∼ 0 to 1, respectively (neglecting variance in Edd or dependence on SMBH mass). For our model, the per-AGN EMRI rates for the 10 6 M SMBH vary from emri ∼ 10 −7 − 10 −5 yr −1 for an SMBH accretion rate Edd = 0.1. Using the above, we derive a ( ∼ 0) volumetric rate of = emri SMBHs / BH ∼ 0.5 − 10 yr −1 Gpc −3 , though we caution that the true rate will be more dependent on the fraction of AGN that are near-Eddington. This simple estimate suggests that the in-situ EMRI rate may be comparable to or greater than the rate derived via AGN-disc assisted binary capture.
Note that these rates do not necessarily correspond to the detectable rate, given that the more massive SMBH in our disc models produces EMRIs with GW frequencies too low to be detectable by LISA. Detectable EMRIs via dynamical processes (e.g. two body encounters) in dry nuclei are expected to occur at a rate 10 − 10 3 yr −1 according to Babak et al. (2017). Overall our results suggests an additional boost to the EMRI rate driven by in-situ formation in sufficiently massive accretion discs.

Final distribution of stars and compact remnants
While most stars are disrupted or accreted by the the central SMBH, a portion of the population can be left over after disc dispersal if the migration timescale is longer than the disc lifetime. This occurs for lower mass stars in discs with higher viscosity parameter = 0.1 (or lower densities) for which the accretion lifetime AGN ∼ 10 7 yr. The maximum mass of these remnant populations scales with the mass of the the central SMBH. For the central SMBH masses of {10 6 , 10 7 , 10 8 }M SMBH, stars with masses below {0.1, 2, 50}M survive the disc phase, respectively. These ranges may change if the stars accrete substantially during their evolution, since an increase in mass will affect the migration rate. Accounting for more complex migration scenarios, such as the possibility of migration traps or, even more likely, the mutual gravitational interactions between AGN stars concurrently with disc-driven migration, would diminish the importance of conditions at birth for the stellar population. Further- Figure 10. EMRI rate per AGN, for discs with Edd = 0.1 with a star formation efficiency = 0.01. Solid symbols show rates for discs with viscosity parameter = 0.1, hollow symbols for discs with = 0.01. We include the EMRI rate for discs models with Type I migration (diamonds), Type I migration plus stellar accretion (hexagons) gap-corrected migration (circles), gap-corrected migration plus stellar accretion (squares). All models assume a star formation efficiency = 0.01, except for the gray circles which adopt = 0.3. more, given that we only consider a single AGN activity cycle, one should be careful not to compare our results directly to observed stellar populations (e.g. stellar discs around Sgr A*) around SMBHs that may have a more complex accretion history.

Tidal disruption
When determining the in-situ EMRI rate, we exclude 'disrupted' stars since the outcome of disruption depends on the stellar structure. It is likely that rather than being promptly disrupted, these stars experience Roche Lobe overflow ) which may or may not be stable, depending on the impact of the accretion disc. For the 10 6 M SMBH discs in particular, surviving tidal forces from the central SMBH is the main obstacle for in-situ stellar populations. These stars might be relevant for the production of Quasi-periodic eruptions (Metzger et al. 2022) or some version of TDEs 13 , which we plan to investigate in a follow-up paper. If these stars survive partial disruptions, they may continue to migrate into the inner disc and potentially increase the in-situ EMRI rate.

DISCUSSION
In this work we show that gravitational instability in radiatively efficient, massive accretion discs produces unique stellar populations. Subsequent evolution of embedded stars can lead to three outcomes: tidal destruction of stars, accretion of stars (or their remnants) as  Table 2. System parameters, stellar population characteristics, and EMRI rates for each subsequent disc and migration model. From left to right we show the central SMBH mass, the disc viscosity parameter, the SMBH accretion rate as a fraction of the Eddington rate, the radius of transition t from the inner stable disc to the fragmenting region, the stellar accretion model (either none or limited Bondi), the peak of the stellar mass distribution, the migration model (Type I or with a gap correction), and finally, the resulting in-situ EMRI rate assuming a star formation efficiency = 0.01 or = 0.3. We disregard models with stellar accretion for high , as this would lead to a large population of stars depleting the gas reservoir. extreme mass ratio inspirals, or leftover stars after disc dispesral. The outcome primarily depends on the interplay of the stellar lifetime and migration timescales, which are dependent on the disc structure. Of course, not all AGN are the same-accretion rates and geometry can vary, and our fiducial model is one of many. We demonstrate that within a single disc model, observables such as the in-situ EMRI rate are dependent on the host SMBH and accretion rate. Our disc model is appropriate for massive discs with near-Eddington inflow rates, or quasars, and not low-luminosity, radiatively inefficient accretion flows which are more common in the local Universe. A proper volumetric rate calculation should take into account the quasar population as a function of SMBH mass and redshift.
Our work complements previous studies on in-situ stellar populations in AGN accretion discs which also suggest that star formation in AGN discs should produce stellar populations with large masses (Goodman & Tan 2004;Levin 2007). One critical difference in our work is that we neglect the interaction between in-situ stars, which may merge with each other. We hypothesize that the realistic outcome is likely a combination of the two. The interaction between stars could reduce the EMRI rate, yet lead to more massive EMRIs. Determining at what rate stars merge or evolve independently will require a thorough investigation of the mutual interactions within the population.
In another sense our model is conservative, given that we adopt a relatively low star formation efficiency for the environment considered ( = 1%) as well as corrections to the initial fragment masses. Both of these ingredients act to reduce the mass of the stellar population. Even under these constraints, our models predict that each accretion episode can produce a combined stellar mass of ∼ 10 2 − 10 4 M during a single isolated accretion episode. SMBHs may experience several accretion episodes, or active duty cycles, during their cosmological evolution (Hopkins et al. 2016), which will have implications for the cosmological EMRI rate.
Several recent works consider how stellar evolution may vary when embedded in a denser, hotter environment than the typical ISM Dittmann et al. 2021;Jermyn et al. 2021Jermyn et al. , 2022. These models predict a range of evolutionary outcomes depending on the initial distance of the star to the central SMBH, the accretion mechanism, and chemical mixing. Their choice of initially solar mass stars is consistent with the range of initial stellar masses derived in this work for discs around 10 6 M SMBHs, however the density in our disc models is considerably higher (by at least 3 or more orders of magnitude), which is beyond the capability of the numerical stellar evolution codes employed in those works. How stellar evolution is altered in higher density environments (along with aspherical accretion or radiative feedback) remains to be explored. The expected outcome based on the models by Cantiello et al. is that stars will accrete rapidly but also eject stellar winds. This ultimately limits their growth to ∼ 10s − 100s M . The constant supply of fuel can lead to a population of 'immortal stars', particularly in the inner disc region . We expect that incorporating this aspect into this work will affect the population of stars that survive tidal disruption by increasing their lifetime, though in reality this will depend on nuances of accretion and feedback. Overall these works highlight the importance of understanding the interplay between stellar evolution, accretion, feedback, and migration.
A strong prediction of gravitational instability is an initial increase in mass inflow due to the generation of structures that increase angular momentum transport. This is supported at different scales: One example at stellar scales arises from observations of FU Ori stars, which are young stellar objects that demonstrate rapid increases in brightness over short timescales (Herbig 1977). These flares are interpreted as boosts in the accretion rate onto the star, which can be explained by their accretion discs becoming gravitationally unstable (Hartmann & Kenyon 1996;Liu et al. 2016). Analogously, we expect that the accretion episodes described in this work should occur with an initially high luminosity that decreases as a portion of the gas is converted into stars. In this case the assumption of a steady-state, constant disc model is not necessarily valid over long timescales, for which the accretion rate will increase once the instability is triggered, then decrease as gas is depleted. Future work will explore how the time dependence of the accretion rate affects the disc geometry and stellar outcomes.
Our work supports the idea that accretion discs in AGN can boost the expected EMRI rate, not only by orbit capture of stars and BHs in the nucleus , but also by in-situ fragmentation. While we consider in-situ EMRIs for a range of SMBH masses, only a subset of these will be detectable by future milliHz detectors. For LISA, the sensitivity lies in the ∼ 10 −4 − 1Hz range (neglecting the details of detector noise as a function of frequency). A circular inspiral generates GWs at a frequency that is twice the orbital frequency, which is at most GW = Ω/ ≈ 4 × 10 −3 Hz ( BH /10 6 M ) 1/2 (3 S / ) 3/2 for the final inspiral around a non-spinning SMBH. For the higher SMBH mass of 10 8 M , a circular EMRI will merge at lower frequencies which fall outside of the mHz band. LISA will be ideal for detecting EMRIs around less massive host SMBHs, for which in-situ generation is critically dependent on whether stars survive disruption. In models where stars grow via a limited Bondi accretion prescription (reaching modest stellar masses ∼ 10M ), this occurs at a rate that is comparable to that expected for EMRIs formed via dynamical estimates.
Gas-embedded inspirals present an exciting class of GW sources for future milliHz detectors. Whether they are captured by the disc via orbit intersections or born in the environment in situ, their environmental evolution leads to distinguishable orbital characteristics: These events will have typically lower eccentricity compared to EMRIs formed via purely dynamical processes and possible spin alignment with the central SMBH (depending on the alignment of both constituent BHs with the accretion plane). If several events are detected, the population properties will be skewed by their accretion and migration history. Individual events also present the exciting possibility of detecting signatures of their environment in the GW signal Yunes et al. 2011;Barausse et al. 2014;Derdzinski et al. 2019Derdzinski et al. , 2021Zwick et al. 2021;Amaro-Seoane et al. 2022), which will allow us to constrain properties of AGN discs in regions that are electromagnetically unresolvable.

Caveats and future considerations
There are many aspects of the picture that deserve further exploration.
Our model assumes a simplified, marginally stable disc where = 0 , which implicitly assumes that some feedback mechanism supports the disc against further collapse. In reality, will have radial and temporal dependence that must be solved consistently with various forms of feedback, e.g. from star formation or accretion onto embedded objects. (see e.g. the feedback-dominated accretion model by Gilbaum & Stone 2021). This also means that the SF efficiency should have a radial dependence (e.g. as in the model of Thompson et al. (2005) or the simulations by Mapelli et al. 2012), which will affect the resulting stellar mass distribution. Such considerations may be critical for understanding the SF in discs around less massive SMBHs, for which fragmentation occurs within the transition from optically thick to thin gas (see the opacity profile of the 10 6 M disc model).
Simulations of fragmenting proto-planetary discs by Boley et al. (2010) and Cha & Nayakshin (2011) find that clumps form on eccentric orbits. If stars are born on initially eccentric orbits, this will affect the subsequent accretion rate (more so if it is dependent on relative velocity with the local gas, as in Bondi) and orbital evolution. We plan to consider in more detail the initial conditions of stellar orbits and their effect on subsequent evolution in future work.
Furthermore, the stellar distribution in the disc that arises in our model will be relevant only for a limited time. As the gas content of the disc diminishes due to fragmentation and accretion onto the fragments and the central AGN, stars and BHs can still evolve as a collection of objects mutually interacting via gravity (Boley et al. 2010;Zhu et al. 2012), essentially as a collisional system. Two-body relaxation and resonant relaxation would occur, as they do in star clusters, and potentially generate more EMRIs. This will have to be investigated in a future work, though, as the phase space structure of the collisional disc that would have to be modelled is not a trivial extension of the case of star clusters, which is widely studied in the literature.
We have considered only one stage of gravitational instability. Over longer timescales, new gas could flow into the nucleus and replenish the outer regions of the AGN disc, triggering another episode of fragmentation. This will be more likely to happen at higher redshift, as more frequent interactions between the host galaxy and smaller galactic companions will trigger gas inflows (Callegari et al. 2011;Capelo et al. 2015). Dynamical instabilities in the nuclear region of the host galaxy, such as nuclear bars, could also trigger repeated episodes of gas inflows (Hopkins & Quataert 2011;. At a minimum, such episodes would occur on a duty cycle of order the orbital time of the nuclear region. The latter timescale is in the range 10 6 − 10 7 yr at scales of 100 pc −1 kpc for a gas-rich galaxy (e.g. Hopkins et al. 2016). Hence our EMRI rate estimate should be taken as a conservative estimate, as several cycles of in-situ star formation such as that modelled here could occur during the lifetime of a massive gas-rich galaxy. It is only between these cycles of gas replenishment that the system will evolve as a collisional gravitational system with a lesser role of gas.
Finally, a word of caution should also be spent about the way migration is treated in our model. Indeed, we are assuming the standard picture of linear torques, but it is known that disc-driven torques are considerably more complicated. For example, stochastic torques arise in both self-gravitating discs (Baruteau et al. 2011) and turbulent, magnetized discs (Nelson 2005). Migration rates can also be affected by the local thermal feedback from the star or an accreting BH (Hankla et al. 2020;Velasco Romero & Masset 2020) or mechanical feedback from stellar winds (Gruzinov et al. 2020;Li et al. 2020). Furthermore, changes to the initial mass distribution, such as those resulting from the inclusion of magnetic fields (e.g. Deng et al. 2020 will affect the subsequent migration rates.

CONCLUSIONS
We present self-consistent models of geometrically thin, steady-state, and near-Eddington accretion discs around SMBHs. By calculating the disc structure and corresponding cooling rates, we determine the conditions for fragmentation and resulting protostar properties. For a set of SMBH masses, disc viscosities, accretion rates, and migration estimates, we quantify the initial stellar mass distribution and a resulting EMRI rate. Our results are summarised as follows: (i) For the parameters considered, discs become gravitationally unstable at sub-parsec distances. Rapid cooling and fragmentation into protostars occurs at a few ×10 −3 − 10 −2 pc from the central SMBH, and extends up to ∼ 0.1pc.
(ii) Initial stellar mass distributions are top-heavy and peak at 0.1 − 1M for discs around BH = 10 6 M , 1 − 10M for BH = 10 7 M , and 10 − 10 2 M for BH = 10 8 M . Assuming a SF efficiency of 1 percent, fragmentation results in * 10 2 − 10 3 stars with a relatively flat mass distribution slope of / ∝ −0.7 * . (iii) Subsequent accretion onto stars via a limited Bondi prescription (which preserves the steady state assumption of the global accretion disc) results in stars increasing their masses by ∼ 2 orders of magnitude over the AGN disc lifetime of ∼ 10 7 yr.
(iv) We calculate the timescale for stars to migrate from their initial positions to the 10 Schwarzschild radii, where they are either tidally disrupted or accreted by the central SMBH. Migration of the stellar population leads to three outcomes: tidal disruption of stars, in-situ extreme mass ratio inspirals with the central SMBH, or leftover stars at sub-parsec separations from the central SMBH.
(v) The rate of in-situ EMRIs produced via the orbital evolution of these stellar populations is emri ∼ 0 − 10 −4 yr −1 per AGN per SMBH accretion episode. Not all models considered produce EMRIs: for BH = 10 6 M , stars need to survive tidal disruption, so EMRIs only occur when stars are able to grow to 10M . Our rate is an under-estimate for SMBHs that experience multiple activity cycles.
(vi) EMRI rates are dependent on the SF efficiency, stellar accretion, and migration prescription. Rates will increase if star formation produces a higher number of massive stars, or if deviations from migration occur that slow down the migration of stars towards radii where they become victim to tidal forces.
(vii) Our disc models suggest the presence of migration traps due to changes in density and temperature gradients. These occur at ∼ 50 − 300 S , depending on the strength of viscosity ( ). At smaller separations, the orbital evolution of stars/BHs in these regions will also be governed by GW emission, and thus the traps only slow down their inevitable inspiral. For the more distant traps, the delay in inward migration may be more substantial.
AGN-born star populations may be diverse depending on the host SMBH and the gas supply, suggesting that this process will vary throughout cosmological history as SMBHs grow.

ACKNOWLEDGEMENTS
We are grateful to the referee, Doug Lin, for his insightful comments that improved this work. We thank Zoltán Haiman, Nicholas Stone, and Alex Dittmann for comments on the manuscript and Tassos Fragos, Lorenz Zwick, Mudit Garg, and Simona Pacuraru for inspiring conversations. The authors acknowledge support from the Swiss National Science Foundation Grant 200020_192092. AD acknowledges support from the Tomalla Foundation for Gravity Research. We also acknowledge use of the following software: NumPy (van der Walt et al. 2011), SciPy (Virtanen et al. 2020), and Matplotlib (Hunter 2007).

DATA AVAILABILITY STATEMENT
The calculations and data used in this study are available upon request from AD. Figure A1. The absolute value of the Type I torque profiles normalized by 0 . Red highlighted regions show where the torque is positive (outward) due to changes in the density and temperature gradients. The total torque (solid black line) is a combination of adiabatic and isothermal, depending on the disc cooling rate. [Top panel:] For M6 and = 0.01: In the inner regions where gas is optically thick, radiative timescales are long and migration is in the adiabatic regime, resulting in a migration trap at ∼ 300 . [Bottom panel]: M6 and = 0.1. the disc density is lower, so the total torque is in between the adiabatic and isothermal regimes. The migration trap occurs at smaller radii (∼ 50 S ) due to changes in the cooling rate. The location of migration traps depends on . Figure A2. Absolute value of the total, isothermal, and adiabatic torque profiles for BH = 10 8 M and BH = 10 7 M for discs with = 0.1. Migration traps occur at ∼ 130 S (for 10 8 M ) and at ∼ 80 S (for 10 7 M ).