Mergers of double NSs with one high-spin component: brighter kilonovae and fallback accretion, weaker gravitational waves

Neutron star (NS) mergers where both stars have negligible spins are commonly considered as the most likely “standard” case. In globular clusters, however, the majority of NSs have been spun up to millisecond (ms) periods and, based on observed systems, we estimate that a non-negligible fraction of all double NS mergers ( ∼ 4 ± 2 %) contains one component with a spin of a (few) ms. We use the Lagrangian numerical relativity code SPHINCS BSSN to simulate mergers where one star has no spin and the other has a dimensionless spin parameter of χ = 0 . 5. Such mergers exhibit several distinct signatures compared to irrotational cases. They form only one, very pronounced spiral arm and they dynamically eject an order of magnitude more mass of unshocked material at the original, very low electron fraction. One can therefore expect particularly bright, red kilonovae. Overall, the spinning case collisions are substantially less violent and they eject smaller amounts of shock-generated semi-relativistic material. Therefore, the ejecta produce a weaker blue/UV kilonova precursor signal, but — since the total amount is larger — brighter kilonova afterglows months after the merger. The spinning cases also have significantly more fallback accretion and thus could power late-time X-ray flares. Since the post-merger remnant loses energy and angular momentum significantly less efficiently to gravitational waves, such systems can delay a potential collapse to a black hole and are therefore candidates for merger-triggered gamma-ray bursts with longer emission time scales.


INTRODUCTION
Nearly exactly 100 years after their prediction (Einstein 1916), gravitational waves (GWs) were detected for the first time in 2015: the LIGO detectors recorded the "chirp" signal from a merging binary black hole (Abbott et al. 2016).The first NS merger event, GW170817, followed soon after (Abbott et al. 2017b) and in its aftermath an electromagnetic firework was observed all across the spectrum (Abbott et al. 2017c).This first gravitational wave-based multi-messenger detection brought major leaps forward for many long-standing problems.For example, it made clear that NS mergers can indeed produce (short) gamma-ray bursts (sGRBs) (Goldstein et al. 2017;Savchenko et al. 2017), the delay between the gravitational wave peak and the sGRB was used to very tightly ⋆ E-mail: stephan.rosswog@astro.su.se constrain the propagation speed of GWs to the speed of light (Abbott et al. 2017b) and the event was also used for an independent estimate of the Hubble parameter (Abbott et al. 2017a).It further helped to constrain the tidal deformability of NSs and thereby some nuclear matter equations of state could be ruled out (Abbott et al. 2017c).Last, but not least, it revealed the ultimate final destiny of a massive binary star system (van den Heuvel & De Loore 1973;Tauris & van den Heuvel 2023), and it confirmed (Tanvir et al. 2017;Abbott et al. 2017c;Cowperthwaite et al. 2017;Smartt et al. 2017;Kasen et al. 2017;Rosswog et al. 2018) the long-held suspicion that NS mergers are major sources of r-process elements in the cosmos (Lattimer & Schramm 1974;Symbalisty & Schramm 1982;Eichler et al. 1989;Rosswog et al. 1998;Rosswog et al. 1999;Freiburghaus et al. 1999).This first event splendidly demonstrated the power and promise of multi-messenger astrophysics for the near future.At the end of the LVC science run O3, ∼ 90 gravitational wave events have been observed (Abbott et al. 2021) and this number will increase rapidly with future science runs at improved sensitivity.While the initially most conservative expectations for NS mergers were component masses near 1.4 M ⊙ , nearly exactly circular orbits (Peters & Mathews 1963;Peters & Mathews 1964) and irrotational NS binaries (Bildsten & Cutler 1992;Kochanek 1992), observations during approximately the last decade revealed a much broader diversity.In terms of masses, NSs with masses down to 1.17 M ⊙ were observed (Martinez et al. 2015) while the highest secure NS masses exceed 2 M ⊙ (Antoniadis et al. 2013;Fonseca et al. 2021).Even more extreme objects that could potentially (but don't have to) be NSs have been detected by the LIGO/Virgo collaboration (GW190814 and GW200210 092254; Abbott et al. (2020Abbott et al. ( , 2023))).In both cases a ∼ 23 M ⊙ black hole (BH) merged with a ∼ 2.6 M ⊙ object.The latter objects could themselves be the result of a binary NS merger (Bartos et al. 2023a).Most recently, another peculiar binary system, PSR J0514-4002E, has been observed that is of particular interest for our study here (Barr et al. 2024): it contains a millisecond pulsar (period ∼ 5.6 ms) in an eccentric orbit (e = 0.71) around a massive compact object with a mass between 2.1 and 2.7 M ⊙ , yet another candidate for a NSs with extremely large mass.Apart from its very large overall mass, this binary system is remarkable, because it could be an example of a double NS (DNS) binary that contains a millisecond pulsar (MSP), i.e. it could be an example of the systems whose mergers we are studying here.This binary system is located in the globular cluster NGC 1851 and, in general, more than half of the NSs observed in globular clusters are spinning with periods of a few milliseconds1 .As we will discuss in Sect.5.1, mergers that involve a MSP are not just a merely academic possibility, but observations suggest that they should account for a noticeable fraction (∼ 4%) of the merging NS binaries.Despite the preference for irrotational systems, the effect of NS spins on the merger process has been explored since the late 1990s, at that time in Newtonian-plus-GW-backreaction simulations.For example Ruffert et al. (1996); Zhuge et al. (1996); Ruffert et al. (1997); Rosswog et al. (1999) explored, in addition to irrotational binaries, also corotating binaries and cases with both NS spins being anti-alligned with the orbital angular momentum.In Rosswog et al. (2000) binaries systems with only one spinning star were explored, though not with a clear motivation based on stellar evolution grounds.More recently, the effects of spins were also explored in full numerical relativity simulations.For example, Kastaun et al. (2013) used approximate initial conditions to study the prompt collapse to a BH and its resulting spin while Bernuzzi et al. (2014) started their general relativistic simulations from constraint satisfying initial conditions mergers, each time both NSs had the same spin.Motivated by NSs in globular clusters, East et al. (2016) modeled dynamical capture binaries with non-zero eccentricities and spins and thereby also considered cases with only one spinning star.More recently, several studies have further explored NS mergers with non-zero spins, see for example Dietrich et al. (2017); Ruiz et al. (2019); Tsokaros et al. (2019); East et al. (2019); Most et al. (2019Most et al. ( , 2021)); Papenfort et al. (2022) and Dudi et al. (2022).Motivated by the insight from binary stellar evolution that mergers that contain a MSP occur at a non-negligible rate, see Sect.5.1, and by the observation of PSR J0514-4002E (Barr et al. 2024), we study here how mergers with one highly spinning component differ in various potentially observable signatures from irrotational sys-tems of the same mass.Our paper is organized as follows.In Sec. 2 we discuss how DNS systems with only one rapidly spinning stellar component may form.In Sec. 3 we summarize our computational methodology, in particular the Lagrangian numerical relativity code SPHINCS BSSN and how we construct constraint-satisfying initial conditions with the code SPHINCS ID that can be linked to either the LORENE (Gourgoulhon et al. 2001;LORENE 2001) or the FUKA (Papenfort et al. 2021; Frankfurt University/Kadath Initial Data solver 2023) initial data solver library.In Sec. 4 we discuss the simulation setup (4.1), the merger morphology (4.2), the gravitational wave emission (4.3), the dynamic ejecta (4.4) and the resulting electromagnetic emission (4.6).Sec. 5 is dedicated to a discussion of our results while the most salient features of NS mergers with a single, rapid spin are summarized in Sec. 6.

FORMATION OF A DOUBLE NS SYSTEM WITH A RAPIDLY SPINNING COMPONENT
The formation of DNS systems, as expected from binary star evolution of a pair of massive stars in isolation, involves a long sequence of stellar interactions (e.g.Tauris et al. 2017;Tauris & van den Heuvel 2023).These interactions include mass transfer and tides, and the binary system must survive two supernova (SN) explosions to remain bound.The NS spins are a crucial testimony of the origin of the progenitor binary.It has been demonstrated that for systems produced in the Galactic disk, the first-formed (and thus recycled) NS in DNS systems can only be spun-up to spin periods longer than ∼ 10 − 12 ms (Tauris et al. 2015) because of the short-lasting and inefficient mass-transfer stages of the progenitor high-mass X-ray binary.
However, the situation is quite different for DNS systems located in dense stellar environments like globular clusters.Here, because of the possibility of close encounters among stars and other binaries, NSs can be assembled in pairs that include a fully recycled MSP via exchange reactions.An example of a DNS candidate system formed via such dynamical interactions is PSR J1807−2500B (Lynch et al. 2012), which is a binary 4.2 ms pulsar found in the globular cluster NGC 6544.
There are currently about 23 DNSs known in our Galaxy, about half of which will merge within a Hubble time (Tauris & van den Heuvel 2023).It is noteworthy that two out of the three globular cluster DNS pulsars, J1807−2500B and J0514−4002A, have been fully recycled to spin periods of 4.2 ms and 5.0 ms, respectively.This is most likely a result of long-term recycling in a low-mass X-ray binary (LMXB) system (Alpar et al. 1982;Radhakrishnan & Srinivasan 1982), later followed by an encounter exchange of the MSP into the present DNS system, see Fig. 1.The bottleneck, in terms of temporal evolution, of the formation channel depicted in Fig. 1 is the nuclear fusion timescale of the low-mass star that evolves to become the donor star in the LMXB system in which the NS is recycled to ms spin period.This donor star should have a mass of less than 2 M ⊙ to secure efficient recycling (Tauris & Savonije 1999;Tauris et al. 2012).The nuclear evolution timescale is thus at least 1.3 Gyr.However, this timescale is often significantly less than the Hubble time and thus plays no big role for the detection statistics of DNS mergers with a fast-spinning NS component given that the in-spiral timescale of the DNS system due to GWs can take any value between a few Myr and up to a Hubble time.
A recent discovery in the Galactic globular cluster NGC 1851 (Barr et al. 2024) provides further evidence for exchange interac- tions producing a MSP orbited by another NS or BH, and may thus hint at the possibility for producing GW source components in the BH lower-mass gap.These components themselves may originate from DNS mergers in dense clusters (Barr et al. 2024) or in hierarchical triple systems (Ransom et al. 2014).To summarize, it is evident that nature may produce DNS systems (and mergers) in which one component is a MSP with very rapid spin.We note that both globular clusters, nuclear star clusters, and hierarchical triple systems have been suggested in relation to the origin of some recently detected GW source components (Rodriguez et al. 2016;Liu & Lai 2021;Bartos et al. 2023b).
The fastest spinning MSP known to date has a spin period of 1.40 ms (PSR J1748−2446ad, Hessels et al. 2006).Recycling an MSP to a spin period of e.g.1.0 ms, requires accretion of, at least, 0.24 M ⊙ of material (see Eq. ( 14) in Tauris et al. 2012), and it is uncertain whether disk-magnetospheric conditions in the vicinity of the accreting NS in an LMXB allows for such fast spin (see discussions in Tauris & van den Heuvel 2023).Nevertheless, in the following we study the merger process of a DNS system in which one NS component has a spin parameter χ ≡ S/M 2 ADM = 0.5, which, somewhat depending on the EOS, corresponds to a spin period of ∼ 1.2 ms (see Table 1), while the other component has a negligible spin (in our simulations χ 2 = 0).The latter approximates the typi-cal spin period ≳ 1.0 s of an old non-recycled NS which has spun down due to significant magneto-dipole radiation from an (initially) strong B-field of order ∼ 10 13 G.
There can be delay times of several Gyr from the formation of the DNS system until it merges, but this is strongly depending on the orbital period (P orb ) and eccentricity of the system, and also depending on the chirp mass (M ); see Peters (1964).For a binary with relatively small eccentricity e ≲ 0.6, this delay time is approximately given by (Maggiore 2008): 1 − e 2 7/2 . (1) As an example, a circular DNS system with an orbital period of 8 hr and two NS component masses of 1.3 M ⊙ will merge after a time interval of about 2.0 Gyr, but the same system with eccentricity e = 0.9 would merge within about 6 Myr.For the circular case, the non-recycled NS would have spun down to a very slow spin period (justifying the use of χ = 0), whereas a recycled MSP with P ≃ 1.2 ms and a low B-field (∼10 7 G) has such a small value of Ṗ ∼ 10 −21 that its position would remain "frozen" in the (P, Ṗ)-diagram and thus it would retain its rapid spin for several Gyr, therefore justifying the use of χ = 0.5.

METHODOLOGY
Here we concisely summarize our simulation methodology as implemented in the current "version 1." of our fully general relativistic, Lagrangian hydrodynamics code SPHINCS BSSN (Rosswog & Diener 2021;Diener et al. 2022;Rosswog et al. 2022Rosswog et al. , 2023)).This novel numerical relativity code evolves the spacetime by solving the full set of Einstein equations via a standard BSSN-formulation (Shibata & Nakamura 1995;Baumgarte & Shapiro 1999).The spacetime evolution part is performed with finite difference methods and with fixed mesh-refinement in a very similar way as in Eulerian numerical relativity codes.Differently from them, we evolve matter with freely moving particles according to a modern formulation of relativistic Smooth Particle Hydrodynamics (SPH).Our approach with a mesh for the spacetime and particles for the matter evolution requires a continuous back-and-forth mapping from the particles to the grid (the mapped property is the energy-momentum tensor T µν ) and from the grid to the particle positions (mapped are the metric properties).Thus, our evolution code consists of three main sectors: matter evolution, spacetime evolution and the coupling between the two.These different elements will be briefly summarized below.

Hydrodynamics
We solve the relativistic hydrodynamics equations with a modern version of Smooth Particle Hydrodynamics, see Monaghan (2005); Rosswog (2009); Price (2012); Rosswog (2015b) for general reviews of the SPH method.A step-by-step derivation of the general relativistic SPH equations can be found in Sec.4.2 of Rosswog (2009), a concise summary is provided in Rosswog (2010) and Rosswog (2015b) and the exact version of the equations that we use together with all technical details can be found in our recent code paper devoted to "version 1.0" of SPHINCS BSSN (Rosswog et al. 2023).
We distinguish between a "computing frame" in which the simulations are performed and a local fluid rest frame in which, for example, thermodynamic quantities are measured.Our equations can be derived from a discretized relativistic Lagrangian2 where ν b is the (conserved) baryon number of SPH particle b, is a generalized Lorentz factor with v being the coordinate velocity, u is the internal energy per baryon and n and s are the local rest-frame baryon number density and specific entropy.We measure a computing frame baryon number density, N = √ −gΘn, g being the determinant of the spacetime metric, in a similar way to Newtonian SPH where the smoothing length h a characterizes the support size of the SPH smoothing kernel W , see below.The equations we evolve in time are the canonical energy per baryon (see Rosswog (2009) for details), and the canonical momentum per baryon where E = 1 + u + P/n is the relativistic enthalpy per baryon with P being the gas pressure.The corresponding evolution equations are where We have implemented a large variety of different kernel functions, but our preferred ones are a C6-smooth Wendland kernel (Wendland 1995) for which we use exactly 300 contributing neighbours, which we had scrutinized in Rosswog (2015a), and a member of the harmoniclike kernels (Cabezon et al. 2008), W H n , with n = 8 for which we use 220 contributing neighbour particles.For the explicit expressions and the motivations behind these choices, we refer to our recent detailed code paper (Rosswog et al. 2023).The hydrodynamic terms are enhanced by dissipative terms to allow for a robust treatment of shocks.See, for example, Fig. 4 in Rosswog & Diener (2021) for a relativistic 3D shock tube test and Sec.2.1.1-2.1.3 in Rosswog et al. (2022) for the explicit expressions for the dissipative terms that we are using.In experiments with the Newtonian MAGMA2 code (Rosswog 2020b) we had found that slope-limited reconstruction in artificial dissipation terms massively suppresses unwanted dissipation effects.We therefore apply this technique also in SPHINCS BSSN.In addition, we also evolve the dissipation parameters in time as described in detail in Sec.2.1.3 in Rosswog et al. (2022).

Spacetime evolution
We evolve the spacetime according to the ("Φ-version" of the) BSSN equations (Shibata & Nakamura 1995;Baumgarte & Shapiro 1999).We have written wrappers around code extracted from the McLachlan thorn of the Einstein Toolkit Einstein Toolkit web page (2020); Löffler et al. (2012), see Sec. 2.3 of Rosswog et al. (2023) for the explicit expressions.The derivatives on the right-hand side of the BSSN equations are evaluated via standard Finite Differencing techniques where we use sixth order differencing as a default.We have recently implemented a fixed mesh refinement for the spacetime evolution, which is described in detail in Diener et al. (2022), to which we refer the interested reader.For the gauge choices, we use a variant of "1+log-slicing", together with a variant of the "Γ-driver" shift condition (Alcubierre et al. 2003;Alcubierre 2008;Baumgarte & Shapiro 2010).

Coupling between the particles and the mesh
The SPHINCS BSSN approach of evolving the spacetime on a mesh and the matter fluid via particles requires a continuous information exchange between the two entities: the gravity part of the particle evolution is driven by derivatives of the metric, see Eqs. ( 6) and ( 7), which are known on the mesh, while the energy momentum tensor that is needed as a source in the BSSN equations is known on the particles.This bi-directional information exchange is needed at every Runge-Kutta substep or, in our case with an optimal 3 rd order Runge-Kutta algorithm (Gottlieb & Shu 1998), three times per numerical time step.The mesh-to-particle step is performed via 5 th order Hermite interpolation, that we have developed (Rosswog & Diener 2021) in extension of the work of Timmes & Swesty (2000).Contrary to a standard Lagrange polynomial interpolation, the Hermite interpolation guarantees that the metric remains twice differentiable as particles pass from one grid cell to another and thus avoids the introduction of additional noise.Our approach is explained in detail in Sec.2.4 of Rosswog & Diener (2021) to which we refer the interested reader.For the more challenging particle-to mesh mapping, we use in the latest version of SPHINCS BSSN a mapping technique that is based on a "local regression estimate" (LRE).This new approach is explained in detail in Sec.2.4 of Rosswog et al. (2023).Here we only summarize the main ideas and refer to our original paper for details.The major idea is to assume that the function to be mapped is known at the particle positions and that this function can be expanded in a Taylor series around any given grid point.This Taylor expansion can be interpreted as an expansion in a polynomial basis with to-be-determined "optimal" coefficients.These coefficients are obtained by minimizing an error functional similar to common least square approaches.In principle, such expansions could be performed up to arbitrarily high order, but (i) the size of the matrices to be inverted and therefore the computational effort increases rapidly with polynomial order, (ii) with higher order the matrices can get closer to being singular and (iii) high order is not everywhere the best approach.For example, when encountering a sharp transition such as the NS surface, high-order polynomials can lead to spurious, "Gibbs-phenomenon"-like oscillations and in such a situation a lower polynomial order delivers a better result.Therefore, we apply a "Multi-dimensional Optimal Order Detection" (MOOD) approach.The idea is to perform the LRE-mapping for different polynomial orders (in practice we use up to quartic polynomials which requires a 35 × 35 matrix to be inverted), discard those results that are considered "unphysical" (e.g.negative energy density, T 00 < 0), and choose out of the remaining options the one which is the best representation of the surrounding particles according to an error measure.This method works very well in practice and we refer the interested reader to Sec. 2.4 and Appendix D of Rosswog et al. (2023) for more details and tests.

Constructing initial conditions
We set up our initial SPH configurations with our initial data code SPHINCS ID, see Rosswog et al. (2023), Sec. 3, for the description of the latest version.This code can be linked to either the LORENE (Gourgoulhon et al. 2001;LORENE 2001) or the FUKA (Papenfort et al. 2021; Frankfurt University/Kadath Initial Data solver 2023) initial data solver library.These initial data solvers provide constraint satisfying solutions for the matter distribution and the corresponding spacetime for a given physical system.This is the first step towards initial conditions for SPHINCS BSSN, but more needs to be done: since we want to use equal mass particles (for purely numerical reasons), the particle distribution needs to accurately reflect the matter density solution (found by the initial data solver) and the distribution of the particles should also ensure a good SPH interpolation accuracy.To achieve this, we use a general relativistic version of the "Artificial Pressure Method" (APM) that was originally originally developed in a Newtonian context (Rosswog 2020a).The main idea is to start from a guess distribution of particles, measure their local error compared to the result found by the initial data libraries and translate this error into an "artificial pressure".This latter pressure is then used in an equation very similar to a hydrodynamic momentum equation to push the particle into a position where it minimizes its error.Our earlier versions of the APM (Rosswog & Diener 2021;Diener et al. 2022;Rosswog et al. 2022) were straight forward translations of the Newtonian method where the "artificial pressure" is calculated from the density error.Recently (Rosswog et al. 2023), we realized that we get -for the same computational effort-slightly more acurate results if we use the error in the physical pressure (rather than in the density) to calculate the "artificial pressure".For details of this latest version we refer to Sec. 3 in Rosswog et al. (2023).While LORENE has been a major work horse for many groups in setting up initial conditions, it struggles in achieving more extreme mass ratios and it does not allow for the construction of binaries with arbitrary spins.In this study, we therefore use initial conditions exclusively based on the FUKA library.

Equations of state (EOS)
Our hydrodynamic equations as described in Sec.3.1 still need to be closed by an equation of state.Currently, we are using piecewise polytropic approximations to cold nuclear matter equations of state (Read et al. 2009), that are enhanced with an ideal gas-type thermal contribution to both pressure and specific internal energy, a common approach in numerical relativity simulations.For explicit expressions please see Appendix A of Rosswog et al. (2022).To date, we have implemented 14 piecewise polytropic equations of state, but for our purposes here, we restrict ourselves to • SLy (Douchin & Haensel 2001): maximum TOV mass For the tidal deformabilities we have quoted the numbers from Table 1 of Pacilio et al. (2022).For all cases we use the piecewise polytropic fit according to Table III in Read et al. (2009) and we use a thermal polytropic exponent γ th = 1.75 as a default.Given the observed mass of 2.08 +0.07 −0.07 M ⊙ for J0740+6620 (Cromartie et al. 2020) the SLy EOS is still above the 2σ lower bound of 1.94 M ⊙ , but probably too soft and we consider it as a limiting case on the soft side.Concerning the currently "best educated guess" of the maximum NS mass, a number of indirect arguments point to values of ∼ 2.2−2.4M ⊙ (Fryer et al. 2015;Margalit & Metzger 2017;Bauswein et al. 2017;Shibata et al. 2017;Rezzolla et al. 2018;Sarin et al. 2020;?), and a recent Bayesian study (Biswas & Datta 2021) suggests a maximum TOV mass of 2.52 +0.33  −0.29 M ⊙ (see also introduction in Godzieba et al. 2021, and references therein), close to the values of e.g. the APR3.From our investigated EOSs, we therefore consider APR3 as the most realistic one which is also consistent with the findings of Pacilio et al. (2022); Biswas (2022).The MS1b EOS with its very high maximum mass of 2.78 M ⊙ and its large tidal deformability is disfavored by the observation of GW170817 (Abbott et al. 2017b), we therefore consider it as a limiting case on the stiff side.

Simulation setup
In this exploratory study, we restrict ourselves to equal mass binary systems with 1.3 M ⊙ for each star and use (in order of increasing stiffness) the SLy, the APR3 and the MS1b EOS.For each case, we run a reference simulation of an irrotational binary and one where one of the stars has a spin parameter χ 2 = 0.5 which, depending on the EOS, corresponds to spin periods between 1.15 and 1.66 ms, close to the 1.40 ms of PSR J1748−2446ad (Hessels et al. 2006).Both spinning and non-spinning cases are set up by using the FUKA initial data solver (Papenfort et al. 2021).The simulations start (t = 0) from an initial separation of 45 km and are performed with 2 million SPH-particles, and initially 7 mesh refinement levels out to ≈ 2268 km in each coordinate direction.The initial minimum grid spacing is ∆x = 369 m, but as laid out in detail in Rosswog et al. (2023), new refinement levels are added dynamically, when a time step criterion is met.These simulations are summarized in Table 1.To illustrate how large the effect of the spin χ 2 = 0.5 is on the stellar structure, we show in Fig. 2 the particle positions of the initial configuration of our reference run APR3 sspin projected on both the XY and XZ plane.The spinning star is significantly rotationally flattened with the extent in the Z-direction being only about 0.8 of the extent in the XY -plane.

Morphology
The spin has a significant impact on the morphological evolution.Since the spinning star is substantially enlarged in the orbital plane, it is more vulnerable to the tidal field from the more compact non-spinning companion and becomes disrupted into a large tidal Table 1.Summary of the performed simulations.The masses are always 2 × 1.3 M ⊙ , the dimensionless spin parameter of star one is always χ 1 = 0 and we always use 2 million SPH particles, an initial separation of 45 km, 7 initial grid refinement levels with a grid spacing ∆x = 369 m on the finest refinement level.The columns are, in order, 1) the name of the run, 2) the equation of state, 3) spin parameter χ 2 of star 2, 4) angular velocity, 5) spin period, 6) total initial ADM mass, 7) total initial ADM angular momentum, 8) the total radiated energy, 9) the total radiated energy in percentage of the initial ADM mass, 10) the total radiated angular momentum and 11) the total radiated angular momentum in percentage of the initial ADM angular momentum.The units are given in the column headings.tail, see Fig. 3, so that the evolution resembles one of an unequal mass binary system.The rapidly rotating high-density core keeps punching into the forming accretion torus, thereby shock-heating the torus, see the fourth columns (t= 11.4 ms) in Fig. 3.A volume rendering of simulation APR3 sspin is shown in Fig. 4. As can be seen from panels two to four, the shocks driven into the forming torus heat up the debris which, impeded by matter in the orbital plane, is expanding vertically and cause a puffing up of the torus.We show in Fig. 5 the evolution of the peak densities (solid lines) and the minimum value of the lapse function (dashed lines), the non-spinning case is shown in blue, the spinning one in orange.The maximum density for the softest EOS (SLy) is about twice as large as for the stiffest EOS (MS1b).During the inspiral, the peak density and lapse are practically identical between spinning and non-spinning binaries for all cases.The two quantities show a clear anti-correlation in the sense that the lapse reaches a minimum at maximum compression and vice versa.For all EOSs, the collision is substantially more violent in the non-spinning case, where larger densities and lower lapse values are reached and the post-merger oscillations are of larger amplitude and persist for longer in the non-spinning case.As expected, these effects are more pronounced for softer EOSs.

Gravitational wave emission
We have used the WeylScal4 thorn of the Einstein Toolkit (Löffler et al. 2012) to postprocess our simulation data and then analyzed the resulting waveforms using Kuibit (Bozzola 2021).In Fig. 6 we show the dominant ℓ = 2, m = 2 mode for all EOSs.In all cases, we see a small initial transient (the first, slightly too large peak) which is due to imperfect initial conditions and subsequently the expected "chirp signal".The waveforms are aligned in time so that the peak of the amplitude coincides at t = 0. Directly thereafter, the merged object is strongly compressed, with the bulk of matter being close to axial symmetry and therefore the GW-amplitude at this stage is minimal.On bouncing back, a bar-like structure forms and the amplitude increases again sharply.For all EOSs, the non-spinning case is shown in black while the spinning case is shown in orange.The gravitational wave amplitudes are largest for the softest EOS and smallest for the stiffest EOS.In all cases the spinning binary has a weaker signal just after merger than the non-spinning binary.This is consistent with the observation that the collision is generally more violent in the nonspinning case (as seen in the density/lapse evolution, Fig. 5, and also in the ejecta properties, see below).Note, however, that in the case of SLy and MS1b, the GW amplitudes become similar about 5 ms after merger.Also observe that in some cases the amplitude decays, but then starts to increase again.This happens at about 7 ms for APR3 and at 15 ms for MS1b.In Fig. 7, we show the radiated energy (solid curves) and angular momentum (dashed curves) in percentage of the initial ADM values.For all EOSs the nonspinning cases are shown in black while the spinning case is shown in orange.Consistent with the observation that the wave amplitude decreases for stiffer EOS, the overall radiated energy and angular momentum is largest for SLy and smallest for MS1b and consistent with the observation that the non-spinning collision is more violent, the radiated energy and angular momentum is significantly smaller for the spinning case.It is also worth noting that much more energy and angular momentum are radiated after the merger than during the simulated inspiral.The SLy EOS cases also still radiate significantly at the end of the simulations, while in the other cases the loss of angular momentum and energy has largely calmed down.For massive enough binaries, the less efficient GW-emission for spinning cases will increase the lifetime of a central remnant before it collapses to a BH (East et al. 2019;Papenfort et al. 2022).In Fig. 8 we show the spectra of the full waveforms for all the EOSs.For the softer EOSs (SLy and APR3), not only are the peak amplitudes smaller but it also shifts moderately (∆ f < 200 Hz) to lower frequencies for the spinning case compared with the nonspinning case.This is consistent with the findings of Bernuzzi et al. (2014); Dietrich & Ujevic (2017) and East et al. (2019).The shift in frequency is largest for the softest EOS and, in contrast, essentially zero for the stiffest EOS.This is consistent with the density evolution shown Fig. 5: for the softer EOSs the density difference between spinning and non-spinning case is much more pronounced than for stiffer EOSs and the larger peak densities result in higher post-merger GW-frequencies.It is also worth noting that the secondary peaks, visible in the spectra of the non-spinning cases, are suppressed in the spinning case for all EOSs.Observationally, such differences in the post-merger phase are likely too small to be measurable at current sensitivities.However, they may be observable with future gravitational-wave instruments (Abbott et al. 2017;Punturo et al. 2010;Ackley et al. 2020).We note that the shift in peak frequency, could also be an important consideration for studies that link the tidal deformability to the post-merger frequency to infer the presence of a hadron-quark phase transition (e.g., Bauswein et al. 2019).

Dynamic ejecta
Since we only perform simulations for ∼ 20 ms and neither include magnetic fields nor neutrino transport, we focus here on the dynamic ejecta.As a criterion to identify them, we use the "Bernoulli criterion", see e.g.Rezzolla & Zanotti (2013), where U 0 is the time component of a particle's four-velocity and E is the specific enthalpy.To avoid falsely identifying hot matter near the centre as unbound, we apply the Bernoulli criterion only to matter outside of a coordinate radius of 100 (≈ 150 km) together with the additional condition of the radial velocity being positive.In a previous study (Rosswog et al. 2022) we had found good agreement between the Bernoulli and the "geodesic criterion", −U 0 > 1.

Masses
In Fig. 9 we show the evolution of the mass identified as unbound by the above described criterion.Clearly, the NS spin has a major effect and increases the dynamic ejecta masses by roughly an order of magnitude.While the irrotational cases yield practically the same ejecta masses for all EOSs, for the cases with spin there is a significant difference between ejecta masses correlated with the EOS stiffness.This is because the bulk of the ejecta is launched via tidal torques that are much larger for less compact NSs.As a corollary this implies that the bulk of the dynamical ejecta is very close to the original NS β-equilibrium value and will be dominated by electron fractions of ∼ 0.1.This will yield heavy r-process (A > 130) which, in turn, will power a red kilonova component.The ejecta masses and velocities are summarized in Table 2.We have additionally labelled the ejecta as "polar", defined as being within 30 • of the rotation axis and as "equatorial" for the rest.Our spinning cases typically eject an order of magnitude more mass than the irrotational cases.In nature, of course, also the mass ratio

MS1b
Figure 6.The ℓ = 2, m = 2 modes of the gravitational wave strain extracted via the Weyl scalar Ψ 4 .The left panel is for the SLy equation of state, the middle panel is for APR3 and the right panel is for MS1b.In all cases the waveforms are plotted with black lines for the non-spinning case and with orange lines for the spinning case.For ease of comparison the waveforms have been shifted in time so that the peak amplitude is at t = 0.
will have an impact on the ejecta masses.Papenfort et al. (2022) interestingly find the largest amounts of ejecta for equal mass cases with large spins (rather than, say, for extreme mass ratios).

Velocities
In Fig. 10 we show the masses in different velocity bins where the velocities refer to those asymptotically reached at infinity.The mass-averaged ejecta velocities hardly exceed 0.2 c, but in each of the simulations we find ejecta velocities exceeding 0.5 c.While this material is not very well resolved, there seems to be a robust trend of the non-spinning cases producing more high-velocity material.This is a result of the less violent collision in the spinning cases discussed in Sec.4.2 where the central object avoids particularly deep compressions which, when bouncing back, produce high-velocity ejecta.

Secular ejecta
Apart from the dynamic ejecta, mergers also produce secular ejecta that become unbound on time scales much longer than what we simulate here.They amount to a substantial fraction (∼ 0.3 -0.4) of the initial disk masses and expand typically at ∼ 0.1 c (Siegel & Metzger 2017;Miller et al. 2019;Fernandez et al. 2019).To estimate the disk masses, we monitor the mass with ρ < 10 13 gcm −3 that is still bound as a function of time, see Fig. 11.Here we also see a substantial increase of the disk masses with spin, however, much less pronounced than for the dynamic ejecta.We find a factor of ∼ two for the softest EOS (SLy) and a substantially smaller increase for the medium (APR3; 25 %) and the stiffest EOS (MS1b; 22 %).The estimated amount of secular ejecta, conservatively estimated as 0.3 M disk , is shown in the last column of Table 2.It is worth stressing that a) the dynamical and secular ejecta masses of the more realistic EOSs (SLy and APR3) add up for the irrotational cases to numbers that are close to but slightly larger than the esti- dE, χ=0 dJ, χ=0 dE, χ=0.5 dJ, χ=0.5

MS1b
Figure 7.The radiated energy (left axis, solid lines) and angular momentum (right axis, dashed lines) in percentage of the initial ADM values.The nonspinning cases are plotted with black lines and the spinning cases with orange lines.The left panel is for the SLy equation of state, the middle panel is for APR3 and the right panel is for MS1b.For ease of comparison the waveforms have been shifted in time so that the peak amplitude is at t = 0.In all case that black line is for the non-spinning case and the orange line is for the spinning case.The positions of the main peaks of the spectra after merger are indicated with vertical dashed lines.mates of ∼0.04 M ⊙ from GW170817 (Cowperthwaite et al. 2017).However, for the spinning cases they can be up to a factor of 2 larger, and b) the increase is the more pronounced for more softer EOS, see Table 2.A substantial disk mass increase through initial spins has also been reported by Papenfort et al. (2022).

Kilonova
With these large differences in ejecta masses due to the spin, one can also expect significant differences in the bolometric light curves of the resulting kilonovae.We will begin by calculating the kilonova emission for only the dynamic ejecta, and, in a second step, we add the secular ejecta whose properties we estimate based on the disk masses.But before we do so, we want to start with simple order of magnitude estimates (Arnett 1980(Arnett , 1982; Kasen & Bild-Table 2. Dynamic ejecta: columns 2 to 7; masses are in units of 10 −3 M ⊙ , the numbers in brackets indicate which percentage is in "polar"/"equatorial" ejecta.All velocities given in units of c. Secular ejecta: columns 8 and 9; column shows 8 the disk masses at the end of each simulation, these numbers should be considered as lower limits.The last column shows an estimate of the secular ejecta msses (= 0.  sten 2010; Metzger et al. 2010;Kasen et al. 2013;Piran et al. 2013;Grossman et al. 2014) to understand the qualitative dependence on various parameters.Let us assume that the ejecta are spherical and have already reached the homologous expansion stage (this happens quickly, see e.g.Rosswog et al. (2014); Neuweiler et al. (2023)) and the radius evolves as R = vt, where v is the expansion velocity and t the expansion time.For an ejecta mass m ej we then have a density of The typical optical depth is τ ∼ κρR, where κ is a typical opacity (keep in mind that in reality this is of course wavelengthdependent).The diffusion time is then given by The peak emission occurs when the expansion time, t = R/v, equals approximately the diffusion time τ diff .After solving for the time, one finds Any initial thermal energy has at τ peak been lost to expansion and thereafter the energy source is radioactive decay where η is some efficiency (neutrinos, for example, are lost) and q is the nuclear heating rate per unit mass.The nuclear heating rate can often be well parametrized as a powerlaw (Metzger et al. 2010;Korobkin et al. 2012;Lippuner & Roberts 2015;Hotokezaka et al. 2017;Rosswog et al. 2018;Hotokezaka & Nakar 2020 with α ∼ 1.3, so that the luminosity is and, at peak, Dynamic ejecta.Since our current set of simulations does not include neutrino transport (or any approximation to it), we have no simulation-based information of the neutron-richness of the ejecta.However, at least for the spinning cases, the ejecta are heavily dominated by tidal ejecta that are never substantially heated and therefore are ejected with their original, cold β-equilibrium Y e < 0.1 (see, for example, Fig. 21 in Farouqi et al. (2021)).We take an angle of |Θ| < 30 • to approximately divide the "polar" and "equatorial" ejecta and assign them values of Y pol e = 0.3 and Y equ e = 0.1.Based on these properties, we compute kilonova light curves with a semi-analytic eigenmode expansion model (Wollaeger et al. 2018;Rosswog et al. 2018) which is based on Pinto & Eastman (2000).We use opacities selected according to the electron fraction (following Tanaka et al. (2020), Table 1).For the heating due to radioactivity, we use a fit formula (Rosswog & Korobkin 2022) that yields the heating rate as a function of velocity, electron fraction and time.This fit formula is based on a heating rate library that has been produced with the Winnet nuclear reaction network (Winteler 2012) that in turn is based on the BasNet network (Thielemann et al. 2011).The network contains 5831 isotopes from the valley of β-stability to the neutron drip line, starting with nucleons and reaching up to Z = 111.We applied a thermal efficiency that is based on Barnes et al. (2016) with parameters suggested in Metzger (2019).The resulting light curves are shown in Fig. 12.The lightcurve peaks for the spinning cases are reached about a factor of 2 later (at ∼ 6 hours post merger) than for the irrotational cases, broadly consistent with Eq. ( 12).As expected from Eqs. ( 15) and ( 16), the irrotational binaries produce substantially dimmer kilonovae with luminosities at 5 days being about an order of magnitude lower than for the single spinning star case.Among those,  the APR3 EOS case achieves the brightest peak luminosities, likely due to the largest amount of polar ejecta.In the spinning cases, the kilonova lightcurves reach their peak several hours after the nonspinning cases.This finding is broadly consistent with those of Papenfort et al. (2022).In Fig. 14, we show the absolute magnitude lightcurves for the kilonova (only from the dynamical ejecta) for two filters, ztfg and ztfr, computed assuming a blackbody spectral energy distribution for the APR3 equation of state and the same Eigenmode expansion model as described above through REDBACK (Sarin et al. 2023).These broadband lightcurves again illustrate the significant differences in the kilonova brightness and evolution for the spinning case.Dynamic plus secular ejecta.We are not modelling the secular ejecta directly, but since this channel likely produces the largest ejecta fraction, we try to estimate its effects on the resulting emission.To this end we (conservatively) assume that 30% of the disk mass becomes unbound, see last column in Table 2, and escapes at 0.1 c.The electron fraction of these ejecta is not fully settled yet, see Rosswog & Korobkin (2022) Sec.2.2.3 for a discussion and pointers to the original literature, therefore we study one case with Y e below (0.20) and one case above (0.30) the critical electron fraction value of Y crit e ≈ 0.25 (Korobkin et al. 2012;Lippuner & Roberts 2015).The corresponding lightcurves are shown in Fig. 13.Clearly, the impact of the electron fraction (which mostly determines the opacity κ) of the secular ejecta is significant for the peak luminosity (see Eq. ( 16)), the peak time (see Eq. ( 12)) and for the slope of the lightcurve decay.To fix ideas, let us concentrate on our best guess EOS, APR3, shown in the middle panel of Fig. 13, similar statements hold for the other EOSs.Changing the secular electron fraction from 0.2 to 0.3 increases the peak luminosity by a factor of ∼ 3, and also delays the time to peak by about the same factor (compare, for example, the solid blue line with the dashed blue line).The impact of the spin is qualitatively similar (compare e.g. the solid blue with the solid orange line): it increases the peak luminosity by ∼ 40 % while delaying the time to peak by nearly a factor of 2. Since the secular ejecta dominate the masses by at least factors of a few, and in some cases more than an order of magnitude, their properties also predominantly shape the electromagnetic emission, so that the effects of dynamical ejecta alone will be very difficult to infer from observations.To again allow a more direct comparison to observations, we also show the absolute magnitude for the combination of secular and dynamical ejecta in Fig. 15.This follows the same trend as the dynamical ejecta only lightcurve, in so far as the spinning case produces a brighter kilonova.However as with the bolometric lightcurve the properties of the secular ejecta shape the electromagnetic emission.

Kilonova precursor and afterglow
We also briefly mention the fast velocity ejecta component and its observational consequences.In each of the simulated cases, we find dynamic ejecta velocities above 0.5 c.Such a fast component may cause additional electromagnetic signatures: in the leading, fastest ejecta, neutrons may avoid being captured by nuclei and their subsequent decay will trigger a blue/ultra-violet precursor to the kilonova ∼ 1 h after the merger (Metzger et al. 2015).Additionally, these fast ejecta, once being decelerated by the interstellar medium, may cause a "kilonova afterglow" due to synchrotron  emission (Nakar & Piran 2011;Hotokezaka et al. 2015Hotokezaka et al. , 2018)).This may have actually been observed as a late-time increase in X-ray flux years after GW170817 (Hajela et al. 2022;Troja et al. 2020;Hajela et al. 2022).To explore the effect of spin on this "kilonova afterglow" in detail, we calculate the kilonova afterglow signature from the dynamical ejecta using two complementary approaches to check for consistency.In the first approach, we take the kinetic energy and total mass for the spinning and non-spinning case numerical simulations for the APR3 equation of state, and model the synchrotron emission following Sarin et al. (2022).In particular, we model the dynamical ejecta as a one-zone spherical outflow propagating out into a constant-density interstellar medium with the observed emission a product of synchrotron radiation produced by the ejecta and external medium (Nakar & Piran 2011) corrected for synchrotron-self absorption and evolution into the deep-Newtonian regime (Sironi & Giannios 2013).As an alternative, second approach, we take the distribution of mass and velocity from the numerical simulations described above, approximating the mass vs velocity relationship as broken power laws and calculate the synchrotron emission following Sadeh et al. (2023).Our simulated kilonova afterglow for a fiducial choice of afterglow parameters; interstellar medium density, n = 5 cm −3 , electron power law index, p = 2.1, magnetic and electron energy fractions of ε e = 0.1, and ε B = 0.001, respectively, at a luminosity distance of 40 Mpc are shown in Fig. 16.
For both modelling assumptions we see the same trend, that the zero spin system peaks significantly earlier than the brighter spinning case.This is consistent with physical intuition, given the latter (spinning) case has more ejected material and a larger kinetic energy.This suggests that kilonova afterglows may provide a potential late-time distinguishing feature for NSs with a spinning component, something we discuss in more detail in Sec. 5. Our numerical modelling agrees with physical intuition, peaking at ≈ 600 and 1200 days for the non-spinning and spinning cases respectively, broadly consistent with the deceleration timescale for our chosen parameters, i.e., the timescale where the blastwave has swept up a comparable mass to its own (Nakar & Piran 2011).

Fallback accretion
While the question whether fallback accretion can power some emission on time scales substantially longer than the dynamical time scale (∼ ms) is clearly very relevant for understanding the ob-     served late time emission of some GRBs, it actually is very difficult to deliver reliable estimates from full-fledged numerical relativity simulations that end at several 10 ms.We will therefore apply different methods to ensure that our conclusions are qualitatively robust.
First, we aim at estimating the mass that will fall back.This is straight forward to calculate in a very simple picture where some mass is launched, follows a ballistic trajectory and will later return to a radius where its available energy is dissipated.The situation at the end of a numerical relativity simulation, however, is more complicated in the sense that the mass distribution is rather continuous so that idealized objects like "the disk" are not straight forward to identify.To illustrate the mass distributions in our runs, we show in Fig. 17 the (fractional) baryonic mass that is included in a given radius.Not too surprisingly, the non-spinning cases are more compact and reach a large fraction of the total enclosed mass, say, 99% already at ∼ 100 km, while this takes substantially larger radii (∼ 600 km) for the χ = 0.5-cases.To estimate the fallback mass, we only consider fluid parts with densities below a threshold of ρ thresh = 10 9 g cm −3 .This value has been chosen after careful inspection of all the simulations (at t = 22.9 ms).Details may depend on the exact value chosen, but none of the main results depends on the exact value.
To ensure the robustness of the estimate, we use two different criteria.First, we use the general relativistic Bernoulli criterion, see Eq. ( 9).Since the specific enthalpy, E = 1 + u + P/n, tends to unity at infinity and the zero-component of the four-velocity tends to the negative Lorentz factor, −Γ ∞ , the condition B ≡ −EU t > 1 means that a matter portion still has a finite velocity at infinity, i.e. it is unbound 3 .If instead B < 1, the flow portion will return and it is considered as fallback.As a second criterion, we use the Newtonian eccentricity of a particle a (e.g.Shapiro & Teukolsky (1983)) where E a , J a and µ a are the particle's Newtonian orbital energy, angular momentum and reduced mass and M the enclosed mass.For this criterion, we identify as "fallback" matter with ρ < ρ thresh and with e < 1.The fallback masses found for both criteria agree very well, even in the worst case the difference is only a few percent, see Table 3.The amount of fallback matter is typically a factor of two larger for the spinning cases.Such fallback could release energies E fb ∼ 2 × 10 51 erg ε 0.1 M fb 0.01 M ⊙ .We also estimate the fallback luminosities based on a simple analytical model (Rosswog 2007).To this end we assume a particle follows a Keplerian orbit from its current position to apocentre and then back to its "circularisation radius", R circ , see e.g.Frank et al. (2002).This is the radius of a circular orbit corresponding to a specific angular momentum value and it is an estimate of the size of a forming accretion disk.We denote the time to reach this radius as t fb , it can be calculated analytically, see Rosswog (2007).After falling back to R circ , it will still take approximately a viscous accretion disk time scale to be accreted onto the central object.This time scale is given by (Shakura & Sunyaev 1973;Frank et al. 2002) Here α is the Shakura-Sunyaev dissipation parameter, ω K , the Keplerian angular velocity and H the disk scale height.So the time from the current position to being accreted is is estimated as t fb +t visc .A similar approach had been taken before in Rosswog et al. (2013).
For the plot of the mass fallback rate, Fig. 18, we use α = 0.01 and H = R circ /3.The spinning cases provide a substantially larger fallback luminosity than the non-spinning cases.They in particular show a pronounced peak in the fallback rate at t ∼ 1 s, where they are about an order of magnitude brighter than the non-spinning cases.The spinning cases also show a longer fallback time scale, at least in this simple model.

DISCUSSION
5.1 How likely are DNS mergers with a rapidly spinning NS?
Unlike the situation for double BH binaries, the fraction of DNS binaries (and mixed BH+NS binaries) formed dynamically in globular clusters has been suggested to be quite insignificant.Ye et al. (2020) estimate the merger-rate density in the local Universe to be ∼ 0.02 Gpc −3 yr −1 for both DNS and BH+NS binaries in globular clusters, or a total of ∼ 0.04 Gpc −3 yr −1 for both populations.In comparison, a conservative merger-rate density estimate based on simulated Galactic field (isolated) DNS and BH+NS systems combined is between 50 and 200 Gpc −3 yr −1 (Kruckow et al. 2018), i.e. four orders of magnitude larger.There are, however, a number of circumstances pointing to a much higher fraction of DNS mergers with a rapidly spinning NS component, including: i) the fraction of DNS mergers with a rapidly spinning NS in the Milky Way is of order 0.04 (see below); and ii) triple systems may also contribute to the formation rate of rapidly spinning NSs in DNS systems (Hamers & Thompson 2019).
Evaluating the fraction of DNS mergers with a high-spin NS component is non-trivial for several reasons.Besides the observed data, one has to factor in the lifetimes of MSPs versus less, or mildly, recycled NSs that are usually observed in DNS systems, beaming morphology and radio survey statistics (Lorimer & Kramer 2012).Note also that dynamical exchange encounter events have a dependence on e.g.orbital separation and cluster age.Thus, a detailed analysis is much beyond the scope of this paper.We will, nevertheless, try here to provide a rough estimate based on a few simple arguments.Four out of the total of 23 DNS systems discovered so far4 in our Galaxy are located within a globular cluster.Three of these systems (PSR J0514−4002A, PSR J0514−4002E and PSR J1807-2459B), corresponding to ∼ 13% of the full DNS sample, host a MSP, whereas the fourth system (B2127+11C), hosting a mildly recycled pulsar spinning at 30 ms, is the only one of those four systems that merges within a Hubble time (it will merge in 217 Myr).Fully recycled MSPs (required for a high-spin NS component investigated here) have, on average, weaker B-fields and are therefore known to have longer spin-down timescales than the mildly recycled NSs that are usually detected in DNS systems.Hence, to correct for this difference in detectability, the current set of empirical data of DNS systems containing an MSP has to be corrected by a smaller production rate (and hence smaller merger rate).From the current sample of DNS systems in the Galactic field, 25% of the mildly recycled NSs in DNS systems have spin-down timescales of τ > 10 Gyr, whereas 30% of the 99 MSPs with spin periods, P < 6 ms in the Galactic field5 have τ > 10 Gyr.Thus, the difference in lifetimes as detectable radio pulsars is, perhaps, not as pronounced as one would immediately guess.However, the average and median values of τ for the 16 DNS systems in the Galactic field, in which a mildly recycled NS is detected, are 4.4 Gyr and 1.0 Gyr, respectively.For the 99 fully recycled MSPs, the values are 9.4 Gyr and 6.0 Gyr, respectively.These values reduce the fraction of DNS systems containing an MSP to 6.1% and 2.2%, respectively.
Among the full population of Galactic binary MSP systems, there is no correlation between the spin of the recycled pulsar and the orbital period up to about 200 days (Tauris et al. 2012).Assuming the same holds for globular cluster DNS systems, i.e. that the formation of DNS systems that will merge is independent of spin period, then, based on the above small number statistics, we may roughly expect a fraction of 4 ± 2% of all DNS mergers to contain an MSP.We emphasize that a more detailed estimate must include e.g.survey selection effects etc. (Kim et al. 2003;Pol et al. 2019).The construction of the full Square Kilometer Array is anticipated to increase the number of known Galactic DNS systems by almost an order of magnitude (Keane et al. 2015).Thus in the coming decade, empirical measures will constrain much better the ratio of DNS systems with one fast-spinning NS component.

Nucleosynthesis
The decompression of NS matter is very appealing as an r-process scenario, since the extreme neutron-richness inside the original stars allows for an effortless and robust production of heavy rprocess (A > 130) up to and beyond the "platinum peak" at A = 195 (Lattimer et al. 1977;Rosswog et al. 1998;Rosswog et al. 1999;Freiburghaus et al. 1999;Korobkin et al. 2012).In the last decade, however, increasingly more ejection channels were identified where the electron fraction becomes substantially increased compared to the pristine value inside the original NSs, see e.g.Wanajo et al. (2014); Perego et al. (2014); Just et al. (2015); Martin et al. (2015); Wu et al. (2016); Siegel & Metzger (2017, 2018); Miller et al. (2019);Fujibayashi et al. (2020a,b).The raised electron fractions can lead to a broad range of r-process elements and this is consistent with the likely identification of strontium in the kilonova of the first detected NS merger GW170817 (Watson et al. 2019).Extremely long-lived central remnants may actually lead to such large electron fractions that the heavy r-process is underproduced.DNS mergers containing a recycled, milli-second pulsar, however, are strongly dominated by tidal ejecta at the original electron fraction and therefore major sources of lanthanides and actinides, a property that they may share with unequal-mass DNS mergers.Since the central remnants are longer-lived due to the less efficient GW-emission, their neutrino emission and absorption by the surrounding debris may lead to an additional, lighter r-process component.Some metal-poor (and therefore old) stars show supersolar ratios of Th or U to Eu (Roederer et al. 2010;Holmbeck et al. 2019) and one may wonder whether a DNS merger containing a ms-pulsar may be a candidate for producing the heavy r-process early on in the Galactic history.The binary evolution leading to the formation of such a binary, however, will take a long time with the bottleneck being the nuclear evolution timescale of the low-mass star that evolves to become the donor star in the LMXB system, see Fig. 1 in which the NS is recycled to ms spin period.Since the nuclear evolution time scale is at least 1.3 Gyr, we do not consider such mergers as good candidates for the very early enrichment of the Galaxy with r-process elements.

Electromagnetic emission
The presence of a rapidly spinning NS can drastically increase the dynamic ejecta, in the studied cases typically by one order of magnitude.Although we do not model weak interactions in the current study, we expect that the vast majority of the dynamic ejecta will be very low in Y e , since it never was shock-heated, and will thus produce a bright red kilonova.This is because a) the initial rotational flattening leads to large tidal ejecta amounts, b) the less violent collision disfavours shock-heated ejecta and c) it also should lead, at least initially, to lower neutrino luminosities disfavoring neutrino captures (and therefore Y e -changes) by the already escaping, far away ejecta matter.Also the secular ejecta are increased, but to a lower degree: we find a factor of ∼ 2 for the softest EOS (SLy) and smaller values for the stiffer EOSs.This ejecta component predominently determines the kilonova properties so that the dynamic ejecta alone will be difficult to infer.The less violent compression and therefore weaker re-bounce in the spinning case also results in a smaller amount of high velocity ejecta.While the latter are only a small fraction of the overall ejecta that is not well enough resolved to draw firm, quantitative conclusions, it seems safe to state that, if indeed they produce a blue transient ∼ 1 hour post-merger (Metzger et al. 2015), this transient should be substantially weaker for the spinning cases.The combination of a weak, blue precursor together with a particularly bright red kilonova may thus be the tell-tale signature for DNS mergers that contain a fast-spinning NS.The high quantity of the ejecta itself could also be a smoking-gun signature of a DNS merger with a fast-spinning NS.However, it is unclear whether a high-spin origin can be robustly estimated from observations, especially in light of other origins of large ejecta amounts such as unequal mass ratios, although these could in principle be verified through gravitational-wave data in the event of a multi-messenger discovery.The kilonova afterglow could also be used to discern whether a merger had at least one spinning NS.We find that mergers with a spinning NS component will produce a brighter kilonova afterglow that peaks at later times, consistent with physical intuition, given the larger mass and kinetic energies.Whether this difference in brightness/peak timescale can be robustly verified given typical uncertainties in afterglow parameters is an open question.While modelling uncertainties in kilonovae and kilonova afterglows may prevent distinguishing a merger with a spinning component from a merger with two non-spinning components, the combination of a weak blue precursor, but brighter kilonova and kilonova afterglow may provide strong constraints on the spin of the NS.For multi-messenger events, such an independent constraint on the spin of the NS through electromagnetic emission, will help to break the degeneracy between component masses and spin and help provide stronger inferences about the mass of merging NSs.We also find that mergers with a spinning component have longer period of fall-back activity.This may provide a more natural solution to explain the duration of a new class of long GRBs that appear to originate from mergers (Rastinejad et al. 2022;Levan et al. 2023) or power late-time X-ray flares such as the low significance flare from GW170817 (Piro et al. 2019).

SUMMARY
In this paper, we have studied a scenario that by and large has probably not been taken seriously enough: the merger of two NSs where one of the components has a spin close to zero while the other is a spun-up MSP.Such binary systems form regularly in regions of large stellar density such as globular clusters and, based on observed systems, we estimate that such binaries may constitute a non-negligible fraction (∼ 4%) of the merging NS binaries, see Sec. 5.1.The merger of such NS binaries has a number of distinctive signatures, that we briefly summarize below, for more details, we refer to the corresponding sections of this paper.
• Morphologically, the spinning (equal mass) cases ressemble mergers with a mass ratio that is significantly different from unity, see Figs. 3 and 4.This is because the spinning stars are substantially more extended due to the rotational flattening, therefore they are more vulnerable to tidal forces and thus produce a single massive tidal tail.
• This implies that the dynamic ejecta are enhanced by an order of magnitude, see Fig. 9, and consist predominantly of low Y e -material at the original NS electron fraction.Also the secular ejecta are enhanced, depending on the EOS, by up to a factor of 2, see Table 2. Spinning cases could -depending on the currently not well determined secular ejecta composition -potentially eject vary large amounts of lanthanide-/actinide-rich matter.They are, however, not good candidates for the source of "actinide boost stars", since the latter must have formed early in the galactic history, while the formation time of the spinning binary NS mergers are determined by the long nuclear evolution timescale of the low-mass star that evolves to become the donor star in the LMXB system.
• The resulting kilonovae are substantially brighter and peak later than in the non-spinning cases, see Sec. 4.6.1 and Figs. 12, 13 and 14.Not too surprisingly, the overall kilonova appearance is determined by the secular ejecta properties, since they constitute the largest ejecta fraction.
• Due to the larger tides of the rotationally flattened spinning stars, the collisions themselves are less violent, lead to smaller compressions of the remnants and to substantially reduced gravitational wave emission, see Table 1 and Fig. 7.
• For the softer EOSs (SLy and APR3) the polar ejecta are about 10% faster than the equatorial ejecta, while for the very stiff MS1b EOS the polar ejecta are actually slower, probably because it is difficult to have shocks with such a stiff equation of state.These statements apply to both the spinning and non-spinning cases.
• The smaller compression in the spinning cases also goes along with smaller amounts of fast ejecta.Therefore, spinning cases produce weaker blue kilonova precursor emission.However, since the ejecta carry more mass and kinetic energy, the kilonova afterglow is brighter and peaks substantially later than in the corresponding irrotational cases, see Fig. 16.
• Since more matter is launched into non-circular orbits, the fallback accretion is substantially brighter and is expected to last longer, see Fig. 18.Whether such systems could be behind GRBs that appear to have a bright kilonova, but are otherwise uncomfortably long for a NS merger (Rastinejad et al. 2022;Levan et al. 2023), given their shorter accretion timescales (Narayan et al. 1992), needs to be explored in future studies.
• As a corollary of the less efficient GW-emission, see Fig. 7, the merger remnant contains larger amounts of angular momentum, should therefore be longer-lived and could potentially power a GRB for a longer time scale.Since the collision is less violent, neutrino driven winds (at least initially) may be weaker, therefore providing an easier path to the launch of ultra-relativistic outflows.
• Last, but not least, given that the rate of such mergers is nonnegligible, caution should be applied in interpretation of GW-data and BNS merger populations where often the prior assumption of negligible spins is considered as the most likely case.

Figure 1 .
Figure 1.Formation scenario of a DNS merger, with one component being a rapidly spinning MSP, produced in a dense stellar environment via exchange encounters.Acronyms used in this figure: NS: neutron star; LMXB: low-mass X-ray binary; MSP: millisecond pulsar; WD: white dwarf; DNS: double NS; BH: black hole.

Figure 2 .Figure 3 .
Figure 2. Initial configuration of simulation APR3 sspin with χ 2 = 0.5.Shown are the SPH particle positions, projected on the XY (left) and the XZ plane (right).The spinning star (each time the one to the right), is substantially larger in the orbital plane (left panel) and substantially rotationally flattened (right panel).

Figure 5 .
Figure5.Maximum densities and minimum lapse values (dashed lines) for the different EOSs.The initial peak densities are 8.41 × 10 14 g cm −3 for SLy, 7.38 × 10 14 g cm −3 for APR3 and 4.32 × 10 14 g cm −3 for MS1b), the density values in the plot are scaled by a factor of 5 × 10 −16 .The non-spinning cases are each time shown in blue, the spinning cases in orange.

Figure 8 .
Figure 8. spectra of the gravitational wave signal.The left panel is for the SLy equation of state, the middle panel for APR3 and the right panel for MS1b.In all case that black line is for the non-spinning case and the orange line is for the spinning case.The positions of the main peaks of the spectra after merger are indicated with vertical dashed lines.

Figure 9 .
Figure 9. Mass identified as unbound as a function of time (masses in units of M ⊙ ; see main text for more information).

Figure 10 .
Figure 10.Ejecta masses (in M ⊙ ) binned according to the velocities at infinity for three EOSs.Results for the irrotational cases are shown in blue, the spinning cases in orange.

Figure 11 .
Figure11.Evolution of the disk masses (masses in units of M ⊙ ), where we use the mass with ρ < 10 13 gcm −3 that is bound to the central object as an estimate of the disk mass.

Figure 13 .
Figure 13.Bolometric kilonova light curves from all ejecta.i.e. dynamic plus secular ejecta, where the latter masses have been estimated from the disk masses.Since the electron fraction of the secular ejecta is not a settled topic, we explore a case with the secular ejecta Y e being below (Y e = 0.2) and another one with being above (Y e = 0.3) the critical electron fraction value of Y crit e ≈ 0.25.

Figure 15 .
Figure 15.Absolute magnitude lightcurves for the spinning (red) and non-spinning cases (blue) for ztfg (left) and ztfr (right) filters for the dynamical and secular ejecta.

Figure 16 .
Figure 16.Radio (left panel) and X-ray (right panel) synchrotron lightcurves at 3 Ghz and 5 KeV, respectively from the dynamical ejecta for the zero spin (blue curves) and spinning (red curves) for the APR equation of state.The solid and dashed curves represent two different modelling assumptions about the kilonova afterglow, while the gray horizontal shaded band represents the sensitivity limit of the VLA and Chandra telescopes.

Figure 17 .
Figure17.Baryonic mass enclosed inside a given radius (M(r)/M tot at t = 22.9 ms for all considered EOSs.The dashed vertical lines indicate the approximate radii beyond which the eccentricity exceeds a value of unity, i.e. matter is unbound.

Figure 18 .
Figure 18.Mass fallback rates (in units of M ⊙ /s) based on the simple analytical model of Rosswog (2007) for our different simulations.

Table 3 .
Fallback accretion: mass estimates for fallback accretion based on either the general-relativistic Bernoulli (m fb,B ) criterion and based on a Newtonian estimate of the eccentricity (m fb,ecc ), both in units of 10 −3 M ⊙ .