A Theory for Neutron Star and Black Hole Kicks and Induced Spins

Using 20 long-term 3D core-collapse supernova simulations, we find that lower compactness progenitors that explode quasi-spherically due to the short delay to explosion experience smaller neutron star recoil kicks in the ∼100−200 km s−1 range, while higher compactness progenitors that explode later and more aspherically leave neutron stars with kicks in the ∼300−1000 km s−1 range. In addition, we find that these two classes are correlated with the gravitational mass of the neutron star. This correlation suggests that the survival of binary neutron star systems may in part be due to their lower kick speeds. We also find a correlation between the kick and both the mass dipole of the ejecta and the explosion energy. Furthermore, one channel of black hole birth leaves masses of ∼10 M ⊙, is not accompanied by a neutrino-driven explosion, and experiences small kicks. A second channel is through a vigorous explosion that leaves behind a black hole with a mass of ∼3.0 M ⊙ kicked to high speeds. We find that the induced spins of nascent neutron stars range from seconds to ∼10 ms, but do not yet see a significant spin/kick correlation for pulsars. We suggest that if an initial spin biases the explosion direction, a spin/kick correlation would be a common byproduct of the neutrino mechanism of core-collapse supernovae. Finally, the induced spin in explosive black hole formation is likely large and in the collapsar range. This new 3D model suite provides a greatly expanded perspective and appears to explain some observed pulsar properties by default.

These data and inferences require explanation in the context of compact object birth, i.e. in the context of corecollapse supernova explosions (CCSNe).The Blaauw mechanism due to unbinding of progenitor binaries as a result of explosive mass loss (Blaauw 1961;Hills 1970;Hoogerwerf et al. 2001) is quantitatively inadequate, though often subdominantly operative.Models that invoke asymmetrical radio pulsar winds (Harrison & Tademaru 1975) are not compelling (Lai et al. 2001).Models that rely on neutrino jets in the context of rapid rotation at birth (Spruit & Phinney 1998) are not consistent with modern CCSN simulations (Janka 2017;Stockinger et al. 2020;Burrows et al. 2020;Burrows & Vartanyan 2021;Coleman & Burrows 2022).The most comprehensive and well-developed theory for pulsar and black hole kicks at birth involves recoil due to the generically asymmetrical supernova explosion (Burrows & Hayes 1996;Scheck et al. 2006;Burrows et al. 2007;Nordhaus et al. 2010Nordhaus et al. , 2012;;Wongwathanarat et al. 2013;Janka 2017;Gessner & Janka 2018;Müller et al. 2019;Nakamura et al. 2019;Coleman & Burrows 2022) and to the accompanying aspherical neutrino emissions (Fryer & Kusenko 2006;Nagakura et al. 2019b;Rahman et al. 2022), the latter generally being subdominant for neutron star birth but important for black hole birth (Coleman & Burrows 2022).
Asphericities arise in various regions of the supernova progenitor star and the timescales of their effects can vary from a few seconds to hours or days.Simultaneous accretion and explosion leads to interactions between ejecta and accreta, and this phenomenon is the first major physical process that influences the kick velocity experienced by the central object.When accretion terminates and the explosive shock wave has swept away most of the outer envelope a few seconds later, this process ceases and the core has achieved its final asymptotic kick state.This marks the end of kick evolution if the explosion is sufficiently energetic, i.e., if the total energy of the outgoing matter is higher than the binding energy of the outer envelope.For less energetic explosions, part of the ejected matter may slow down due to interactions with stellar envelopes, and may fall back onto the central object (Schrøder et al. 2018;Chan et al. 2020).This second process may last hours to days.Even if the progenitor fails to explode and there is no asymmetric ejecta, the asymmetric neutrino emission can still provide a modest kick (Fryer & Kusenko 2006;Nagakura et al. 2019b;Rahman et al. 2022).In this case, the accretion of the outer envelope can carry angular momenta to the compact core (Antoni & Quataert 2023).
In this paper, we focus on the first process mentioned above, i.e., kicks during the simultaneous accretion and explosion phase.This is the dominant process in energetic CCSN explosions (Coleman & Burrows 2022;Chan et al. 2020).Many relevant simulations have been done in the past (Burrows & Hayes 1996;Scheck et al. 2006;Burrows et al. 2007;Nordhaus et al. 2010Nordhaus et al. , 2012;;Wongwathanarat et al. 2013;Janka 2017;Gessner & Janka 2018;Müller et al. 2019;Nakamura et al. 2019;Coleman & Burrows 2022), but the required sophisticated and expensive 3D simulations have not been continued to late enough times after bounce to achieve asymptotic kick speeds, with very few exceptions (Stockinger et al. 2020;Coleman & Burrows 2022).Most 3D simulations are stopped within one second or less after bounce, far to soon to witness the final kick speed imparted.Moreover, even the exceptional longer-term 3D studies did not explore the systematics with the broad range of progenitor masses to determine the correlations between kick speed and progenitor mass and/or initial core structure, the latter frequently associated with the compactness parameter 1 .They focused on the few models whose recoil kick had quickly settled to nearly a final state.In contrast to papers incorporating detailed numerical simulation, papers such as Janka (2017) and Mandel & Müller (2020) and Bray & Eldridge (2016) and Bray & Eldridge (2018) have derived semi-analytical and/or empirical relations between kick velocities and supernova/progenitor properties.Although such simple recipes are intuitively appealing, they need to be tested thoroughly using 3D numerical simulations.
Hence, a comprehensive theoretical study of the broad spectrum of systematic behaviors in the context of modern 3D core-collapse supernova theory has not yet been performed.With this paper, we attempt to remedy this situation by providing and analyzing twenty of our recent long-term state-of-the-art 3D Fornax (Skinner et al. 2019) simulations for initially non-rotating ZAMS (zero-age main-sequence star) solar-metallicity progenitors from 9 to 60 M ⊙ taken from Sukhbold et al. (2016) and Sukhbold et al. (2018).Our goal is to build a more solid foundation for the current CCSN kick theory and to provide detailed 3D models that can be used to test, calibrate, and update widely-used recipes.
We emphasize that we are addressing in this paper only solar-metallicity progenitors and those most likely to be the origin of most galactic neutron stars and black holes and leave sub-solar metallicity studies to future work.We also are not here addressing the issue of the potential role of very massive massive stars and pulsational pair instabilities2 .Importantly, many of our models have been carried in 3D beyond four seconds after bounce and show signs of asymptoting to their final recoil kick speeds.A few of these models have or will form black holes.
In this paper, we complement our kick study with a study of the associated induced spins (Blondin & Mezzacappa 2007;Rantsiou et al. 2011;Coleman & Burrows 2022).Though a compelling explanation for pulsar and neutron star spins is spin up from an initially spinning Chandrasekhar core upon collapse, modern supernova simulations leave protoneutron star (PNS) cores spinning due to random and stochastic post-explosion accretion of angular momentum (L), even when the core is initially non-rotating.In fact, the magnitude of such induced spins (Coleman & Burrows 2022) can be comparable to measured pulsar spins (Chevalier & Emmering 1986;Lyne & Lorimer 1994;Manchester et al. 2005;Faucher-Giguère & Kaspi 2006;Popov & Turolla 2012;Igoshev & Popov 2013;Noutsos et al. 2013).Moreover, there are hints that radio pulsar spins and kicks might be correlated at birth (Wang et al. 2006(Wang et al. , 2007;;Johnston et al. 2005Johnston et al. , 2007;;Ng & Romani 2007;Noutsos et al. 2013) Finally, angular momentum transport from the stellar core to the mantle up to the time of core collapse might perhaps be more efficient than hitherto thought (Cantiello et al. 2014), leaving the Chandrasekhar core slowly rotating just before implosion.Therefore, we continue here the exploration of the possibility that pulsar spins might be partially reflective of induced spins during the supernova.Similarly, black holes created in supernova explosions (Burrows et al. 2023) will have both large kicks and significant induced spin parameters (a = Lc/GM 2 ), though black holes created in a non-exploding context should have birth spins reflective of the initial progenitor spin (though perhaps in a complicated fashion).However, similar to the case with kicks, the induced spins of the central object can be influenced by different very long-term processes.In this paper, however, we focus on the crucial first few seconds post-bounce.
We used the code Fornax (Skinner et al. 2019;Radice et al. 2019;Nagakura et al. 2019a,b;Vartanyan et al. 2019;Burrows et al. 2019;Burrows et al. 2020;Vartanyan & Burrows 2020;Nagakura et al. 2020;Vartanyan et al. 2022Vartanyan et al. , 2023;;Burrows et al. 2023) to generate the twenty 3D models to late times after bounce used here to study recoil kicks and induced spins in the context both of supernova explosions and quiescent black hole formation.A more comprehensive set of studies using these new 3D simulations is in preparation which will address other aspects and outcomes of core-collapse and their dependence upon progenitor.This is the largest set of long-term (many seconds after bounce) 3D state-of-the-art core-collapse simulations ever created and we hope to use it to gain new insights into the corecollapse supernova phenomenon.To create this set, we employed the CPU machines TACC/Frontera (Stanzione et al. 2020) and ALCF/Theta and the GPU machines ALCF/Polaris and NERSC/Perlmutter.The specific setups and run parameters are described in Burrows et al. (2023) and Vartanyan et al. (2023).The techniques employed to derive the kicks and induced spins are described in detail in Coleman & Burrows (2022) and follow the standard approaches in the literature (Scheck et al. 2006;Wongwathanarat et al. 2013;Stockinger et al. 2020).However, for the calculation of the neutrino component, unlike in the Coleman & Burrows (2022) paper, we calculate its contribution here at the PNS radius (defined as the 10 11 g cm −3 isosurface), and not at the outer mass surface of the final proto-neutron star (PNS).The result for the final total vector kick is the same, but there is a different partitioning between the neutrino and matter components.The momentum transfer between the matter and the radiation in the region between these two definitions accounts for the shift and the total momentum is conserved under both approaches.We employ the different method here to comport with the standard practice in the literature (Stockinger et al. 2020).

KICKS TO NEUTRON STARS AT BIRTH
Until recently, 3D models experiencing kicks had been calculated to post-bounce times insufficient to actually witness their final values.Moreover, due to computational limitations, the systematics with progenitor core structure and total ZAMS mass were obscured, since only one or a few models early during the accumulating recoil could be simulated.This is akin to viewing the entire progenitor panorama through a straw.The stochasticity of the final values of any explosion observable inherent in the chaotic turbulence of the CCSN phenomenon also served to compromise extrapolation and interpretation.However, with our new set of twenty long-term 3D full-physics simulations, we are now able to extract what we suggest are real trends across the progenitor continuum and find correlations among observables hitherto only imperfectly glimpsed.Hence, though for all core-collapse supernova observables there are expected to be a distribution of values, even for the same progenitor, the availability of a raft of such models across the progenitor continuum mitigates in part against this modest impediment to determining system-wide trends.This is the philosophy of our multi-model study, enabled by the efficiency of Fornax and the availability of high-performance computing resources.
Table 1 provides summary numbers of relevance to this study.The top rank displays the data for the exploding models that leave neutron stars (NS), while the bottom rank displays the corresponding data for the models that leave black holes (BH) (see §2.1), either quiescently or explosively.Most of the models have either asymptoted to the final state (in particular the lowest-mass, lowest-compactness progenitors) or nearly so.The left-hand side of Figure 1 depicts the temporal evolution of the absolute value of the vector kick, including both the matter and neutrino recoils and all momentum flux, pressure, and gravitational forces (Coleman & Burrows 2022).The lowest-mass progenitors generally flatten/asymptote earliest, while the models with the largest compactness (frequently more massive) take longer.Such long-duration 3D simulations have not in the past been available.
The right-hand side of Figure 1 depicts the final kicks versus final gravitational mass of the corresponding neutron stars.This is a plot of observables that has not been available before and speaks to the importance of long-term 3D simulations to provide testable predictions for modern CCSN theory.We note that the lowest mass progenitors (blue) are all clustered at the low end of the kick range.This is correlated with the greater degree of sphericity of models that explode early, before the post-shock turbulence can grow to significant strength and Mach number to manifest the growth of large bubbles, that when explosion ensues lead it and set an explosion direction.For these initially non-rotating progenitors, the direction of explosion is random, but once determined a crude axis is set.This is not the case for the lowest compactness progenitors that explode early and, hence, more spherically.
Figure 2 depicts this dichotomy.Its left-hand side portrays the blast as traced by the 10% iso-56 Ni surface of the explosion of the 9 M ⊙ progenitor, colored by Mach number.Purple is positive (outward) and blue is inward (infall).The shock wave is the spherical bluish veil around the 56 Ni bubbles.In contrast, the right-hand side of Figure 2 depicts the corresponding plot for the exploding 23 M ⊙ model3 .This model has been carried to ∼6.2 seconds after bounce and is only slowly asymptoting, while the 9 M ⊙ progenitor flattened within 2 seconds (Coleman & Burrows 2022).The blast structure of the 23 M ⊙ model is more generic, and this is reflected in the separation of the blue dots on the right-hand side of Figure 1 from the rest.The latter are spread between ∼300 km s −1 and ∼1000 km s −1 and the former are found near ∼100−200 km s −1 ; the difference reflects the dichotomy illustrated in Figure 2 4 .Binary neutron stars are typically found to be on average on the low-mass side (Fan et al. 2023) and they seem a distinct population.One explanation is that the lower-mass neutron stars are birthed (on average) with lower kick speeds so that such binary systems can survive.This correlation emerges naturally, without tweaking, from our simulations and is consistent with observations (Tauris et al. 2017).Such a correlation was anticipated by Müller et al. (2019) and explored empirically by Bray & Eldridge (2016) and Bray & Eldridge (2018).The right-hand side of Figure 1 confirms this correlation, and also suggests that there is a neutron star mass threshold, roughly around 1.35 M ⊙ , above which the kick velocities show much larger variations.It is also the case that even lower-mass progenitors not in our sample (such as the 8.8 M ⊙ progenitor of Nomoto (1984) and accretion-induced collapse (AIC) systems) would naturally experience small kicks and these systems could play a role.Such progenitors are also likely to explode quickly and experience only very modest kicks (though if rapidly rotating could experience asymmetric MHD jets that could result in interesting recoils).The lowest progenitor mass range and AIC systems might be relevant as sources for the neutron stars bound in globular clusters, for which the kick speeds need to be smaller that the potential well depths (≤50 km s −1 ).Though we believe the qualitative picture is now clear, we are still some distance from a fully quantitative explanation of the binary neutron star conundrum, as suggested in population synthesis studies.More 3D numerical simulations are required to generate a complete quantitative theory that can fruitfully be compared to the observed mass-kick distribution of neutron stars.Moreover, we emphasize that the progenitor models we use in this work are all single star models (Sukhbold et al. 2016(Sukhbold et al. , 2018)), and that binary effects may significantly influence progenitor structures and, hence, kick velocities (Podsiadlowski et al. 2004).
The left-hand side of Figure 3 depicts the temporal evolution of the corresponding neutrino contribution to the kick.This is its absolute value and in the total kick provided in Figure 1 the matter and neutrino terms are vectorially added; the two are not generally aligned and their individual instantaneous directions vary with time.We see that the lowest ZAMS mass progenitors (such as the 9 M ⊙ model) have some of the highest neutrino kicks (Coleman & Burrows 2022), but that most other models at higher compactness whose explosions are more delayed have lower neutrino kick contributions.Moreover, since the matter contributions to the kicks for the more delayed explosions are generally larger, the neutrino contributions are proportionately less for them and can be as low as ∼3%.They are, on the other hand, as high as ∼80% for the lowest mass/compactness progenitors.We emphasize, however, that the chaos of the CCSN phenomenon can result in a range of values for all observables and that this range has not yet been well determined.
On the right-hand side of Figure 3 we provide a plot of total kick magnitude derived versus the mass dipole of the ejecta.The latter is a measure of the explosion asymmetry.We see that the two are roughly correlated, as expected for the general recoil model of kicks, but that there is a bit of scatter.In Figure 4, we provide a similar plot of the kick magnitude versus the explosion energy.The rough trends of kick with explosion energy and with ejecta mass dipole is clear and generally expected (Burrows et al. 2007;Janka 2017).

Black Hole Kicks
Which progenitor massive stars leave stellar-mass black holes at their endpoints is not yet firm.Despite years of speculation, most of which was uninformed by what at any time was the current (though evolving) state-of-theart supernova theory, there remains no consensus among workers in the field concerning the progenitor/black-hole mapping.The simple idea that above a ZAMS mass cut black holes are the residues of collapse and below that cut neutron stars are left now seems untenable.We ourselves have found (though not proven) that a 60 M ⊙ star leaves a neutron star (see Table 1).This is consonant with the discovery in the young cluster Westerlund 1 of a neutron star (Muno et al. 2006) that from the cluster's turnoff age must be from a progenitor with a ZAMS mass greater than 40 M ⊙ .This also calls into question a helium core mass cut (perhaps near ∼10 M ⊙ ) (Burrows 1987), above which black holes form, despite the fact that the stellar-mass black-hole mass function seems to peak near this mass (MacLeod & Grindlay 2023) and such solar-metallicity stars should often leave stripped-envelope structures (either by mass transfer or winds).
Moreover, a simplistic compactness cut, above which black holes are birthed at stellar death, is also inconsistent with the emerging theory of a turbulence-aided neutrino-driven supernova explosion.Neutrino-driven explosions with energies near ∼10 51 ergs require high compactness to achieve.Without rapid initial rotation and an associated MHDdriven explosion, the low compactness progenitors are not capable of generating explosions in the 10 51 ergs (one Bethe) range.Higher compactness is necessary to provide the higher accretion luminosities and higher "optical" neutrino depths in the absorbing post-shock mantle and gain region (Bethe & Wilson 1985) to provide the power to drive the explosion much above the few ×10 50 ergs we find for the lowest mass progenitors (8.8 to 9.5 M ⊙ ), which possess the lowest compactnesses.
Therefore, a simple mass or compactness cut, above which black holes are the endpoints and below which neutron stars are born doesn't seem correct.This leaves this issue unresolved.However, we (Wang et al. 2022;Tsang et al. 2022;Burrows et al. 2020;Burrows et al. 2023) have recently found a region of progenitor space between ∼12 and 15 M ⊙ where for some subset of initial density profiles in that interval we do not see explosions via the neutrino mechanism using the collection of solar-metallicity massive stars above 9 M ⊙ (Sukhbold et al. 2016(Sukhbold et al. , 2018)).A fraction in this gap skirt explosion in the context of our 2D and 3D simulations, and this fraction could be ∼one-third to ∼two-thirds.This outcome could be because 1) the progenitor models we have employed have some flaw or because 2) our simulations may have some flaw or do not capture some essential physics.Nevertheless, we emphasize that in the context of our detailed 3D simulations the mass interval for this channel is almost wholly a function of the progenitor model suite, which is still very much provisional.Which mass range gives birth to black holes quiescently via the route we identify might shift significantly.
However, perhaps such an interval is indeed a nursery for stellar-mass black hole formation.Two of the models in the second rank of Table 1, the 12.25 and 14 M ⊙ models, are in this interval and do not explode.Even after ∼two to three seconds after bounce, their stalled shocks continue to sink in radius.Curiously, when the shock is this compact, a spiral SASI mode (Blondin & Shaw 2007) is excited that modulates the neutrino and gravitational-wave signatures (Vartanyan et al. 2023) in a diagnostic way.Generally, it is only when the stalled shock sinks below ∼120 km (usually below ∼100 km) that any form of SASI mode emerges from under neutrino-driven convection and this usually presages a failed explosion (though see below).
As Table 1 indicates, the kicks to these black holes are quite small.Since the models don't explode, there would be no matter recoil and the kick would be due solely to anisotropic neutrino emission.If all the mass that resides in these models were to be accreted, given the momentum imparted to the cores and the eventual accretion of the total mass of ∼11−12 M ⊙ remaining after the integrated mass loss up to collapse is factored in (Sukhbold et al. 2016(Sukhbold et al. , 2018)), kicks of ∼7−8 km s −1 are expected.However, we think it quite likely that further mass loss due to possible Roche-lobe overflow, common envelope evolution, and/or impulsive mass loss near the time of collapse (Morozova et al. 2018;Antoni & Quataert 2023) may shave more mass from these progenitors before collapse, leaving ∼8−12 M ⊙ left for black holes created through this channel.Interestingly, this mass range roughly coincides with the peak of the measured range of masses of the relevant X-ray sources (MacLeod & Grindlay 2023).
Note that since these models shouldn't experience the general-relativistic instability until minutes to hours later (Vartanyan et al. 2023), and since we have not simulated to those times, our models manifest some stochastic recoil (see Table 1) due to turbulence around the still accreting PNS.However, we suggest that since the near-terminal stages of massive star evolution and Type IIp light curves (Morozova et al. 2018) show signs of late-time mass ejection that may be aspherical, the recoil momentum from such events might find its way into the black hole.The result could be another component to the speed with which the black hole is born.Whether this component could have a magnitude of ∼100−200 km s −1 (Mandel 2016;Fragos et al. 2009;Repetto et al. 2017;Atri et al. 2019;Mandel et al. 2021;Dashwood Brown et al. 2023) is to be determined.Otherwise, we find that the kicks imparted to such black holes not associated with explosions are below ∼10 km s −1 , with some likely scatter.
However, as we have recently published (Burrows et al. 2023) (see also Chan et al. (2018)), there may be another channel of black hole formation that is accompanied by vigorous explosions.Despite (or because of) the fact that these models have the highest compactness in the Sukhbold et al. ( 2016) solar-metallicity collection, the 19.56 and 40 M ⊙ progenitors exploded vigorously with energies of > 2 × 10 51 ergs and 1.9×10 51 ergs, respectively, as of ∼4 and 21 seconds after bounce.
Curiously, the spiral SASI that emerged for the lower-mass black hole formers that did not encourage explosion did so for their much-higher compactness counterparts.This is because the high compactness models also experience high accretion rates that translate into high neutrino luminosities.This combination of high luminosity at a time that the spiral SASI slightly pushed out the shock and enlarged the gain region was adequate to ignite explosion.As Table 1 and Figure 3 show, this channel of black hole formation not only is explosive with large energies, but has very large kicks.This is due to the anisotropic neutrino-driven jets and anisotropic accretion during explosion, each of which contribute to significant recoil kicks (Burrows et al. 2023).Moreover, for these models, even though the neutrino contribution is still subdominant, it reaches values of hundreds of km s −1 .By the end of these simulations, the recoil momentum achieved is the highest among our set of twenty progenitors.With an estimated final gravitational mass of ∼3.5 M ⊙ (the rest is exploded away), the 40 M ⊙ progenitor is racing at ∼1000 km s −1 .Such a black hole is unlikely to remain bound to a companion.Therefore, we suggest there are at least two channels of black hole formation − one leaving black holes quiescently in the ∼8−11 M ⊙ range with low kick speeds and another leaving black holes explosively in the ∼2.5 − 3.5 M ⊙ mass range with significant proper motions.We note that we have not studied sub-solar-metallicity models, for which the core compactness structures and wind mass loss are quite different, that we have not simulated models above 60 M ⊙ , and that we have not addressed pulsational pair-instability models (Woosley & Heger 2021) likely to provide another, though rarer, black hole formation channel.Rahman et al. (2022) do address that channel.However, they employ a 2D flux-limited diffusion code.Though there are several interesting results in that paper, in our experience (and acknowledged by those authors) 2D codes are compromised in determining kick magnitudes and correlations and 3D models are generally preferred.

INDUCED SPINS
Though the rotation profiles of the cores of massive stars at stellar collapse are not known, we know that massive star surfaces rotate.The loss of angular momentum to winds is guaranteed, but its magnitude is a function of surface B-fields, which themselves are less well-known.Furthermore, the redistribution of angular momentum in the stellar interior is no doubt via magnetic torques, but the processes and couplings are not firmly constrained.The fact that red giant interiors are rotating more slowly (Cantiello et al. 2014) than inferred using a Taylor-Spruit dynamo (Spruit 2002) suggests a more efficient spin-down process than employed by Heger et al. (2005) in estimating the spin at collapse of massive-star models.The latter found using a Taylor-Spruit prescription for angular momentum redistribution that the Chandrasekhar cores might have initial periods at collapse near ∼30−60 seconds, which would translate by angular momentum conservation into neutron star periods of ∼30−60 milliseconds (Ott et al. 2006), with the lowest mass massive stars spinning the slowest.Hence, since these periods are not much different from those expected for radio pulsars (Chevalier & Emmering 1986;Lyne & Lorimer 1994;Manchester et al. 2005;Faucher-Giguère & Kaspi 2006;Popov & Turolla 2012;Igoshev & Popov 2013;Noutsos et al. 2013), it would seem that the issue of the initial spin period of the cores of massive stars is resolved, at least in broad outline.However, it is still possible that many of the cores of massive stars destined to undergo collapse may be torqued down significantly.This thought is motivated not only by the results of Cantiello et al. (2014), but by observation that white dwarfs are born rotating very slowly (Kawaler 2015;Fuller et al. 2019).Could this be the fate of a good fraction of the white dwarf cores of massive stars when they die?In order to explore this possibility, one first needs theoretical insight into what induced spins are possible in the context of core-collapse supernovae.This is what we provide in this paper.
The left-hand side of Figure 5 renders the evolution of the angular momentum imparted to the PNS cores during their post-bounce evolution.This spin up is due to the stochastic accretion of plumes of infalling matter into the PNS, with varying impact parameters.If the explosion is very asymmetrical, as is the case when the kick speeds are significant, or black hole formation is accompanied by an explosion, the consequently asymmetrical accretion during explosion (mostly then from one side) can spin the PNS up to interesting periods.For the black holes that form in an explosion, these spins can reach high values (see Table 1), even though their cores were initially not rotating.This implies that such objects could be collapsar candidates (Heger & Woosley 2003), despite their lack of initial spin and their solar-metallicity progenitors.
The right-hand side of Figure 5 portrays the dependence of the final induced spin period of the neutron star upon ZAMS mass.We see a trend, roughly recapitulated with the kick systematics (see Table 1), of decreasing period with increasing ZAMS mass, but with a lot of scatter.Such scatter should be expected in the context of the chaotic turbulence seen in 3D simulations of CCSNe.We emphasize that the lowest ZAMS mass and compactness progenitors which explode early and approximately spherically have the slowest induced spin rates, while those exploding later and at higher compactness have greater induced spin rates.In fact, our highest compactness models can leave cold neutron stars with high angular momenta and spins.However, while we are comfortable with the compactness/kick and PNS mass/kick correlations we derive, we are much less sanguine that induced spins are the whole story of compact object birth spins.Importantly, we know that the Crab explosion involved a low-mass progenitor near ∼9−10 M ⊙ (Stockinger et al. 2020), but that its birth spin was ∼15 milliseconds (Lyne et al. 2015).This is far from the induced periods (many seconds) we find for this mass range of progenitors, not only suggesting, but requiring, the story to be more complicated.Spin-up by binary interaction is also distinctly possible (Fuller & Lu 2022).Therefore, the relative roles of induced and initial spins in the distribution of observed birth periods for compact objects remains open.
Nevertheless, we find a weak, but intriguing, correlation between induced spin and recoil kicks, demonstrated in Figure 6.The greater the kick, the greater the induced spin, again with significant scatter.Is there something of importance in this apparent correlation?Is there rough support for this, however complicated by the overlapping effect of the initial spin, in the pulsar database?
The induced spins for the non-exploding black hole channel (see Table 1) are quite low, again resulting solely from from the torques due to the radiated neutrinos.Spin parameters of close to zero are inferred from our simulations and model extrapolations.However, should there be any significant angular momentum in the accreting stellar envelope that on long timescales will bulk up the core (something quite likely), these numbers can be radically changed.Hence, even under the assumption that the initial core spin rate is zero, we can not say with any confidence what the final spin parameter of black holes born through this channel, might be.
However, the induced spins for the exploding black hole channel can be large.This can follow from the significant asymmetry in their explosions, itself related to the large recoil kick.A spin parameter as large as ∼0.6 emerges for the 40 M ⊙ model (see Table 1), though we have yet to determine the corresponding number for the 19.56 M ⊙ model.Their large compactnesses and related large post-explosion asymmetrical mass accretion rates are direct causative factors.In addition, accretion continues to power the driving neutrino emissions, which for this black hole formation channel are quite large.Accretion settles into biased accretion onto one general side and spins up the periphery of the PNS before black hole formation.After black hole formation, continued accretion parallels continued spin up until after many seconds the black hole is isolated (for the 40 M ⊙ model after ∼50 seconds?).This general black hole formation scenario is discussed in Burrows et al. (2023).

SPIN-KICK CORRELATION
Figure 2 portrays the Mach number for 56 Ni surfaces and the crudely axial explosions when the compactness is not small and when there is a greater delay to explosion.It also indicates that there is continued accretion roughly in the perpendicular plane.Hence, when the compactnesses and ZAMS masses are not small, explosion and accretion are roughly perpendicular.Accretion then spins up the periphery of the PNS and the net spin vector is either roughly parallel or anti-parallel to the kick (which is along or anti-parallel to the explosion axis).Hence, we would think a spin-kick correlation (Holland-Ashford et al. 2017;Katsuda et al. 2018;Ng & Romani 2007) would naturally emerge.Figure 7 depicts the dot product of the final spin and kick unit vectors for our model set, versus ZAMS mass (left) and compactness (right).We see that though there seems to be a deficit in the perpendicular direction and there is a slight preference for an (anti-)aligned orientation the null hypothesis can not be rejected.Indeed, more rigorously, we perform a Kolmogorov-Smirnov (K-S) test with the null hypothesis that the simulated angles are drawn from a uniform distribution on the sphere.We see a D-statistic of 0.19, which means the cumulative distributions are not very different.In addition, due to the small number of models, the p-value of this K-S test is 0.51, which means that there is a chance that the deficit in the perpendicular direction is simply due to statistical fluctuations.Therefore, however intriguing, many more simulations are required to properly address a potential spin-kick correlation in the initially non-rotating case.
Nevertheless, such a correlation may arise when the Chandrasekhar core is initially rotating, since such rotation may slightly bias (depending upon its magnitude) the axis and ejecta mass dipole vector of the subsequent neutrinodriven explosion.If such is the case, then the spin-kick correlation may emerge.Continued accretion is very roughly perpendicular to the explosion direction, itself roughly (anti-) parallel to the kick direction, resulting in a (anti-)correlation between spin and kick.Note that we have not invoked B-fields nor MHD jets and that this is an expected outcome in the neutrino-driven explosion context.However, we emphasize that we have yet to perform the requisite rotating simulations.

CONCLUSIONS
With a focus on the recoil kicks imparted to the residues of stellar core collapse in the supernova and black hole formation contexts, we have analyzed twenty state-of-the-art 3D long-term core-collapse simulations generated using the code Fornax.In the past, most 3D models were not of sufficient duration to witness the cessation of net acceleration and the asymptoting of the kicks to their final values.This includes our own previous work (Coleman & Burrows 2022), as well as that of other workers in the field.However, and for the first time for a large and uniform collection of 3D supernova models, we have either asymptoted the kicks, or come within 20% of doing so.As a result, we obtain an integrated and wide-angle perspective of the overall dependence of the recoil kicks and induced spins upon progenitor mass and their Chandrasekhar-like core structures, the latter indexed approximately by compactness.As had been postulated in the original recoil kick model (Janka & Mueller 1994;Burrows & Hayes 1996), the recoil is always opposite in direction to the bulk of the ejecta.
We find that lower mass and lower compactness progenitors that explode quasi-spherically due to the short delay to explosion experience on average smaller kicks in the ∼100−200 km s −1 range, while higher mass and higher compactness progenitors that explode later and more aspherically give birth to neutron stars with kicks in the ∼300−1000 km s −1 range.We also find that these two classes can be correlated with the gravitational mass of the residual neutron star.The correlation of a lower kick speed with lower neutron star gravitational mass suggests that the survival of binary neutron star systems, for which their mean mass is lower than average, may in part be due to the lower kick speeds we see for them.With some scatter, there is a correlation of the kick with both the mass dipole of the ejecta and the explosion energy.We also find, as did Coleman & Burrows (2022), that the contribution of the recoil due to anisotropic neutrino emission can be important for the lowest mass/compactness progenitors.
However, unlike Coleman & Burrows (2022), we partition the fractional contribution of the neutrino component to the overall kick differently, which results in an apparent diminution of its contribution for higher ZAMS mass and compactness progenitors, for which the matter recoil effect then predominates.The difference with Coleman & Burrows (2022) is a consequence of our use here of the standard literature radius at which to calculate the neutrino contribution, and not due to some major differences with the overall results in Coleman & Burrows (2022).Between the effective radius at which the neutrino component is calculated in Coleman & Burrows (2022) and here, there is momentum transfer between the matter and neutrino radiation that is properly taken into account in the vector sum of the two provided in both papers.
We find two channels for black hole birth.The first leaves masses of ∼10 M ⊙ , is not accompanied by a neutrinodriven explosion, and experiences small kicks, perhaps below ∼10 km s −1 .The second (Burrows et al. 2023) is through the vigorous, neutrino-driven explosion of a high-compactness progenitor that leaves behind a black hole with a mass of perhaps ∼2.5−3.5 M ⊙ that is kicked to high speeds, perhaps above ∼1000 km s −1 .This second black-hole formation channel is novel and unexpected, but if true challenges most notions of what can explode.Both channels challenge prior notions of the mapping between progenitor and outcome and must still be verified.However, in the context of our CCSN simulations, the associated mass ranges for black hole formation depend entirely upon the as yet unconverged progenitor suite model.In particular, the ZAMS mass range for our quiescent black hole formation channel could shift significantly when binarity, overshoot, shell-merger, mixing processes, and the 12 C(α,γ) rate are properly understood.
Associated with kicks are induced spins that accompany them in the context of the modern theory of turbulenceaided, neutrino-driven supernova explosions.We find that the induced spins of nascent neutron stars range from a few seconds at low compactness and low progenitor mass to ∼10 milliseconds at higher compactness and higher ZAMS mass.However, this can not be the whole story of the birth spins of neutron stars, and the relative roles of initial and induced spins have yet to be determined.Nevertheless, we find that it is in principle possible to explain the some spin rates of some pulsars with the post bounce asymmetrical accretion of net angular momentum in the context of asymmetrical supernova explosions for which spherical symmetry is generically broken.Moreover, we find that a spin/kick correlation for pulsars is not yet suggested by the theory.However, if there is an initial spin to the collapsing core, its vector direction could bias the direction of the explosion.If this is the case, a spin/kick correlation may emerge.However, this speculation has yet to be verified, but, if true, a spin/kick correlation (in a statistical sense) would naturally emerge as a byproduct of the turbulence-aided neutrino mechanism of core-collapse supernova explosions.
In the non-exploding case of black hole formation, the induced spins are quite low and due to anisotropic neutrino emission.Any slight neutrino torques would be due to a random stochastic residual effect.Under these circumstances, an initial spin would be determinative.However, for the exploding context of black hole formation, the induced spin is likely to be large.For our 40 M ⊙ progenitor, we find that the spin parameter could reach ∼0.6 (the maximum is 1.0).This is in the collapsar range, and suggests that even an initially non-rotating, solar-metallicity, high compactness progenitor could be a seat of long-soft gamma-ray bursts, with a vigorous explosive precursor.This is a controversial claim, and one that needs much further scrutiny.However, this is what is suggested by our simulations.
It must be remembered that the chaos in the turbulence generic and fundamental to core collapse phenomena will result in a spread of outcomes, even for a given progenitor.This translates into distribution functions for the explosion energies, nucleosynthesis, morphology, kicks, and induced spins that are not delta functions and are currently unknown.Hence, there is natural scatter in the kicks and induced spins for a given initial model that has yet to be determined, but that will spread observed values.Some of the scatter seen in this study no doubt originates in this behavior.Such are the wages of chaos.Nevertheless, though the study of a single 3D model could distort one's understanding, the study of a wide range of progenitor structures should provide a more complete picture of the theoretical trends and correlations that will minimize some of the confusion due to such scatter.In any case, this is the philosophy of this 3D multi-model, long-term exploration of the core-collapse, but much more remains to be done.In addition, we note that the mapping of observables to progenitor structure, modeled loosely as "compactness," seems more robust than the mapping with ZAMS progenitor mass, tethered as it is to the progenitor suite employed.Given this, the fact that we are in this study using the Sukhbold set of progenitors limits what one can say about the observable/ZAMS-mass correlations and the guidance we can provide to population synthesis modelers.The effects of binarity (mass transfer), overshoot, shell mergers, perturbations, rotation, the 12 C(α,γ) nuclear rate, and wind mass loss have yet to be fully understood and their influence retired.Moreover, physical processes that might change residual black hole properties on timescales of hours, days, or years, such as "fallback" (Schrøder et al. 2018;Chan et al. 2020) or random late-time angular momentum accretion (Antoni & Quataert 2023), are not covered in this work.
Though Fornax is a sophisticated tool, it incorporates no approach to neutrino oscillations, uses a moment-closure approach to the radiative transfer, and incorporates approximate, not precise, general relativity.Nevertheless, the greatly expanded panorama the new capabilities and this new numerical database provide is qualitatively more comprehensive and predictive, due to its breadth and temporal extent, than has heretofore been possible and appears to explain many observed pulsar properties by default.(Breu & Rezzolla 2016).We assume here that the final neutron star radius is RNS = 12km.For black holes, the relativistic spin parameter a = Lc GM 2 is shown and the number in brackets is our guess for its of final value.The final column contains the cosine of the spin/kick relative angle at the end of each simulation.
are the total kick speeds and the kick speeds induced by neutrinos only, respectively.Both are measured at tmax.Numbers in brackets show the estimated final values (kick speed, spin parameter) of the black holes.Note that the 40-M⊙ model has been simulated to 21 seconds after bounce, at which point the spin parameter has reached ∼0.35, and that we estimate its final gravitational mass and spin parameter (a) will be ∼3.5 M⊙ 0.6, respectively, achieved after a further ∼25 seconds of evolution.We have yet to determine the corresponding numbers for the 19.56 M⊙ model.These numbers still need to be checked.L central is the angular momentum of the central object at tmax, and P cold is the inferred final spin period of the cold final neutron star, assuming P = 2πL central I NS .The neutron star moment of inertia is estimated using INS = (0.237 + 0.674ξNS + 4.48ξ 4 NS )M R 2 NS , where ξNS = GM R NS c 2 is the compactness parameter of the neutron star itself

Figure 1 .
Figure 1.Left: Total kick velocity as a function of time in each simulation.Solid lines are models that form neutron stars, while dashed lines are black hole formers.Right: Kick versus cold NS gravitational mass.Black hole formers are not included here.

Figure 2 .
Figure 2.An isosurface (10%) of the 56 Ni abundance in the explosion of a 9 M⊙ (left) and 23 M⊙ (right) progenitor model at ∼1.95 and ∼6.2 seconds, respectively, after bounce, painted by Mach number.Red/purple is positive and white/blue is negative, indicating infall.This figure is a representative depiction of a generic feature of core-collapse explosions by the turbulence-aided neutrino-driven mechanism, mainly that there is often simultaneous explosion in one general direction and accretion in a very roughly perpendicular direction when the explosion time is delayed by a few hundred milliseconds (as for the 23 M⊙ model) and a more spherical explosion when it is launched early (within ∼50−100 milliseocnds), as for an initially non-rotating (or slowly rotating), low-compactness progenitor, such as this 9 M⊙ model.This fact has a bearing on the morphology of the debris, line profiles in supernova ejecta, and possible spin-kick correlations.Since most explosions are delayed and manifest some bubble-driven dipolarity, the explosion direction for most progenitors (anti-)correlates with the kick direction.See text for discussions on this finding.

Figure 3 .
Figure 3. Left: Neutrino kick velocity as a function of time.Solid lines are models that form neutron stars, while dashed lines are black hole formers.Right: Mass dipole versus final total kick speed.Note that those models that seem to deviate from linearity (e.g., the 23 M⊙ model) are also those models which have yet to asymptote, suggesting that a roughly linear relationship is even more robust than this figure implies.

Figure 4 .
Figure4.Magnitude of the final kick speed versus explosion energy.The range of energies reflects the fact that some of these model energies have not yet asymptoted; we provide an informed range when they don't and the simulation of these models is still proceeding.

Figure 5 .
Figure 5. Left: Total angular momentum as a function of time for each simulation.Solid lines are models that form neutron stars, while dashed lines are black hole formers.Right: Estimated neutron star final spin period versus progenitor ZAMS mass calculated assuming L tot I N S .Only neutron star formers are shown here.The neutron star moment of inertia is estimated using (0.237 + 0.674ξ + 4.48ξ 4 )M R 2 NS , where ξ = GM R NS c 2 is the compactness parameter of the neutron star itself(Breu & Rezzolla 2016).We assume here that the final neutron star radius is RNS = 12 km.

Figure 6 .
Figure 6.Relation between kick speeds and induced spins of neutron stars birthed in initially non-rotating progenitor stars.See text for a discussion.

Figure 7 .
Figure 7. Cosine of the angle between the neutron star kick and induced spin directions versus ZAMS mass (left) and compactness (right).

Table 1 .
Results for the simulations presented in this work.MZAMS is the initial ZAMS mass of the progenitor model.tmax is the maximum post-bounce time reached by the simulation.For the black hole formers, the times are the time of black hole formation, and the second time (if shown) is the latest time after black hole formation we have carried that model hydrodynamically.The types of central objects are summarized.NS means neutron star, BH (F) means black hole formed with a failed explosion, and BH (S) means black hole formed accompanied by a successful CCSN explosion.M baryonic is the baryonic mass of the central object at tmax, while Mgrav is the gravitational mass of a cold, catalyzed NS.M final is the estimated final black hole mass.The final masses of the non-explodingblack hole formers are the total masses in the progenitor models at collapse, assuming all matter will eventually be accreted, while the final masses of the exploding black hole formers are derived from our hydro-only blast wave simulations after black hole formation.p kick gives the linear momentum of the object at tmax.v