Evolutionary Channels for the Formation of Double Neutron Stars
We analyze binary population models of double-neutron stars and compare results to the accurately measured orbital periods and eccentricities of the eight known such systems in our Galaxy. In contrast to past similar studies, we especially focus on the dominant evolutionary channels (we identify three); for the first time, we use a detailed understanding of the evolutionary history of three double neutron stars as actual constraints on the population models. We find that the evolutionary constraints derived from the double pulsar are particularly tight, and less than half of the examined models survive the full set of constraints. The top-likelihood surviving models yield constraints on the key binary evolution parameters, but most interestingly reveal (i) the need for electron-capture supernovae from relatively low-mass degenerate, progenitor cores, and (ii) the most likely evolutionary paths for the rest of the known double neutron stars. In particular, we find that J1913+16 likely went through a phase of Case BB mass transfer, and J1906+0746 and J17562251 are consistent with having been formed in electron-capture supernovae.
Since the discovery of the first double neutron star (DNS) system (Hulse & Taylor, 1975) almost 40 years ago, astronomers have discovered only 8 more DNS among more than 2000 pulsars. As exotic endpoints of binary evolution, DNS form through a wide variety of binary phases, including tides, strong stellar winds, two supernovae (SN), multiple phases of mass transfer (stable and unstable), and orbital decay due to gravitational wave radiation (van den Heuvel & Heise, 1972; Flannery & van den Heuvel, 1975; Delgado & Thomas, 1981; Bhattacharya & van den Heuvel, 1991; Tutukov & Yungel’Son, 1993). Binary disruption due to a SN or coalescence during dynamically unstable mass transfer and common-envelope (CE) evolution are both highly likely, yet all observed DNS have survived these phases. Without two supernovae and CE evolution, the observed DNS in their tight orbits could not have formed.
Our current understanding of the full CE evolution and the final outcome for the binaries going through it is incomplete, and consequently CE modeling for binary populations has been limited to relying on simple parameterizations of the problem which involve the amount of mass and energy or angular momentum lost in this non-conservative process (Webbink, 1984; Livio & Soker, 1988; de Kool, 1990; Iben & Livio, 1993; Nelemans et al., 2000; Webbink, 2008; Taam & Ricker, 2010). Although such parameterizations are of limited use in terms of understanding the physics of this phase, they have allowed considerable progress in understanding the formation and evolution of interacting binary populations with compact objects (white dwarfs, neutron stars, black holes). Thus there is always interest in constraining them empirically using observed systems and populations in the Milky Way and other galaxies. Recently, more physically motivated treatments of the CE phase have been put forward (Deloye & Taam, 2010; Ge et al., 2010; Woods & Ivanova, 2011; Ivanova & Chaichenets, 2011; De Marco et al., 2011; Ivanova et al., 2013, for a review). These represent important advancements, however they require the use of a stellar structure and evolution code as part of the modeling, and therefore they are not applicable to the current tools available for population syntheses (which adopt fitting formulae based on single star evolution).
Binaries have a significant chance of getting disrupted when one of the stars goes through a SN, either because of the associated mass loss and/or because of asymmetries in the underlying mechanism imparting a “natal kick” to the newborn NS of hundreds of km/s. Based on observations of proper motions of single pulsars, not only do Hobbs et al. (2005) derive a Maxwellian kick distribution with a dispersion velocity of 265 km/s, but they argue in favor of an absence of pulsars with very low proper motions. However an early population synthesis study by Portegies Zwart & Yungelson (1998) argued that the formation of the observed systems requires a distribution with a substantial fraction of NS born with low velocity kicks. Pfahl et al. (2002b) arrived at the same conclusion, finding that low kick velocities are necessary to explain the observed population of Be/X-ray binaries, and the presence of large numbers of neutron stars in Galactic globular clusters (Pfahl et al., 2002a). In search of a physical explanation for how some neutron stars may form with small natal kicks, Podsiadlowski et al. (2004) drew renewed attention to a core-collapse mechanism known as an electron-capture supernova (ECS). Originally proposed by Miyaji et al. (1980), an ECS is a low energy SN occurring when an ONe degenerate core-collapses due to electron captures onto Mg and Ne (Nomoto, 1984; Miyaji & Nomoto, 1987; Nomoto, 1987). Podsiadlowski et al. (2004) postulated that this channel was only accessible to binary stars that had lost their envelopes in a mass transfer phase, thereby limiting the second dredge up phase. The remnant NS is predicted to receive a kick of km/s. Very soon after the discovery of the double pulsar (Burgay et al., 2003), analysis of its properties led a number of groups (Willems & Kalogera, 2004; Dewi & van den Heuvel, 2004; Piran & Shaviv, 2005; Stairs et al., 2006; Willems et al., 2006; Dewi, 2010; Wong et al., 2010; Ferdman et al., 2013) to argue (albeit for different reasons and at different levels of quantitative analysis) in favor of low-velocity natal kicks being involved in the formation of the second pulsar. In a qualitative analysis van den Heuvel (2007) expanded the picture arguing that through both mechanisms, the standard Fe core-collapse and the ECS, the entire population of DNS could be formed. Most recently Wong et al. (2010) analyzed all individual systems in view of all current observational properties (binary and kinematic) and concluded that some observed DNS must have received high natal kicks and other systems require low natal kicks, also indicating that both mechanisms are needed.
|DNS||ddValues for eccentricity and orbital period are rounded to 3 decimal places.||ddValues for eccentricity and orbital period are rounded to 3 decimal places.||Pulsar Mass||Companion Mass||Evol. ConstraintseeWe only use conservative evolutionary constraints. Most systems do not have enough information to place any constraints on the evolution. See §3.2 and references for a more complete explanation.||References|
|B1913+16||0.617||0.323||59.0||1.4408(3)||1.3873(3)||Channel I or II||2, 3|
|J1906+0746aaThe observed pulsar in J1906+0746 is the unrecycled secondary.||0.085||0.166||144.1||1.248(18)||1.365(18)||9, 10|
|J17532240bbJ17532240 is unconfirmed as a DNS, we therefore exclude it from our analysis. See §4.3||0.304||13.638||95.1||11|
|B2127+11CccB2127+11C is contained within the globular cluster M15 and is suspected to have formed dynamically, therefore we exclude it from our analysis.||0.680||0.335||30.5||1.35(4)||1.36(4)||12|
References. – 1 – Wolszczan (1991), 2 – Hulse & Taylor (1975), 3 – Taylor & Weisberg (1982), 4 – Burgay et al. (2003), 5 – Janssen et al. (2008), 6 – Ferdman et al. (2014), 7 – Corongiu et al. (2004), 8 – Champion et al. (2004), 9 – Lorimer et al. (2006), 10 – Kasian (2008), 11 – Keith et al. (2009), 12 – Jacoby et al. (2006)
In the study presented here we develop population models of DNS formation and evolution and we contrast our results against DNS measurements of orbital period and eccentricity (only). We complement upon the work of previous population synthesis studies of DNS (Kiel et al., 2008, 2010; Osłowski et al., 2011; Dominik et al., 2012) by adding the consideration of any constraints on the evolutionary channel each of the observed systems has followed; we use a Bayesian analysis to quantitatively assess the results. The last analysis to examine whether measured and values are in agreement with DNS population modeling was presented by Dewi (2010), who focused primarily on the formation of J07373039, with no quantitative statistical analysis.
We undertake this study to address the following specific questions: (i) are there current DNS population models that are consistent with DNS and measurements? (ii) can such consistent models also account for the strong evolutionary constraints that are currently available for a few of the observed systems? (iii) if yes, then how do these select models constrain binary evolution processes and parameters; which processes and parameters are not important? (iv) can we use the select models to uncover the evolutionary history of the rest of the observed DNS systems? We choose to restrict ourselves to the and measurements for a number of reasons: they are the best measured DNS properties (without any indirect inference involved), they are highly correlated making it easy to evaluate the models in a 2-dimensional plane, and population models make reliable predictions for these properties, unlike others (e.g., rates). Our additional consideration of evolutionary constraints for specific systems introduces a new element, never used before in DNS population synthesis. In §2 we describe our modeling method followed by our quantitative analysis method in §3. In §4, we present our results and discuss them in §5, ending with a summary in §6.
2. Modeling Methods and Simulations
2.1. StarTrack Population Synthesis Code and Models in This Study
To simulate DNS formation, we adopt the population synthesis code StarTrack (Belczynski et al., 2008). We use a reference model with a given set of parameters, and we then vary a large number of them from the reference assumptions to assess their effect on DNS properties and formation channels.
StarTrack adapts the fitting formulae by Hurley et al. (2000) for single star modeling (time evolution of key macroscopic physical properties) and it modifies them to account for the effects of mass transfer to allow for modeling of binary evolution. The mass-transfer treatment (outside of the case of CE phases) has been tested against and calibrated to mass-transfer sequences (see Belczynski et al., 2008, and references therein). This code has been used for a wide range of studies and applications involving modeling of interacting binaries with all types of compact objects, and has been compared to and demonstrated agreement with observations of such binaries both in our Galaxy and in external galaxies since the early 2000’s when the first version of the code was produced. It has also been maintained and kept up to date as our understanding of compact object formation and stellar winds has evolved over the years. In particular, updated prescriptions for massive stellar winds and their dependence on metallicity have been included (Belczynski et al., 2010). Although a recent revision of the code uses a compact object remnant mass function physically motivated by calculations of the supernova mechanism (Fryer et al., 2012), our results rely on the remnant mass function calculated by Timmes et al. (1996). An additional recent parametrization of the envelope binding energy parameter , a measure of the central concentration of a stellar envelope which is relevant for CE calculations, was added to Startrack after our simulations were completed (Dominik et al., 2012) and is not included in the version of the code used here.
Initial conditions for each population model are as follows: the primary star in the binary is formed on the main sequence with a mass drawn from a power law distribution (Scalo, 1986), while the secondary is determined based on the mass ratio value drawn randomly from a flat distribution between 0 and 1; the initial orbital separation is logarithmically flat and the initial eccentricity follows the dynamical-equilibrium distribution . Allowing for wind mass loss, the primary star is then followed through the Hertzsprung gap (HG), red giant branch, He-main sequence, and asymptotic giant branch (AGB). If the primary is massive enough, it will collapse through a SN, leaving a compact remnant. Assuming the binary is not disrupted during the first SN, the secondary evolves through its lifetime, with possible interruptions by mass transfer sequences and also driven by wind mass loss and/or tides. Once it has evolved completely, if massive enough, the secondary will also leave a compact remnant as the product of a SN. Asymmetries during the SN impart an (assumed instantaneous with respect to the orbital period) kick to the neutron star or black hole, altering the shape of the binary orbit, following the equations of Kalogera (1996, 2000). The formation of a DNS requires that both stars have pre-SN masses within a specific range, such that the binary neither disrupts due to the SN kick nor does it merge prior to the second SN. The particular variables and prescriptions used are given in §2.1.1. The DNS population is given a constant birth rate over the past 10 Gyrs. After the second SN, the orbital evolution due to gravitational radiation is followed using the equations of Peters (1964) until the current epoch.
2.1.1 Reference Model
Our reference model (Model 1) adopts most of the the standard values suggested in Belczynski et al. (2008). The few modifications are described below.
An Fe core-collapse supernova occurs when a star forms a non-degenerate CO core massive enough to ignite stable nuclear burning, creating progressively heavier elements until the degeneracy pressure is overcome and the Fe core burns explosively. If instead the CO core is formed partially degenerate and reaches a critical mass of 1.08 , stable nuclear burning will ensue, creating a degenerate ONe core. We assume that the entire CO core is burned to form the ONe core and approximate that this grows at the same rate that the CO core would. A resulting ONe core more massive than 1.38 is thought to collapse in an electron capture supernova (ECS) while less massive cores form ONe white dwarfs (WD). In our code we vary the range within which the He core mass forms a partially degenerate CO core (). In our reference model, we set this range to 2.0 / 2.5.
To stars undergoing an Fe core-collapse SN, we apply kick velocities randomly chosen from a Maxwellian distribution with a dispersion velocity of 300 km/s, in agreement with observations of the space velocities of single pulsars (Hobbs et al., 2005). However Kitaura et al. (2006) find that the neutrino heating mechanism occurs faster and produces a smaller ejecta mass in an ECS, and they conclude that ECS events impart a smaller kick velocity to collapsing ONe degenerate cores ( km/s). In practice these specific boundary values on the CO core, as well as its progenitor He-rich degenerate core, are quite uncertain. Therefore in our code, we apply a kick decreased by a factor of 10 to a star collapsing in an ECS compared with the standard Fe-core case.
The mass transfer parameters on which we focus in this analysis primarily deal with CE evolution. These phases are treated using the equations provided by Webbink (1984) and de Kool (1990). In this prescription, the efficiency of orbital energy loss during a CE phase is determined by the efficiency parameter, , and the binding energy of an envelope is parametrized with . Since the two parameters are degenerate in the energy conservation formalism for CE evolution, in our models we vary the product, . We set to be 0.5 for the reference model. We do not attempt to calculate the mass accreted within a CE, instead limiting the mass accreted to a random amount between 0.05-0.1 . The fraction of mass lost that is accreted onto a main sequence star during stable Roche lobe overflow is set at 0.5 for the reference model. Stellar evolution codes show that HG stars lack a steep entropy gradient across the core-envelope boundary (Ivanova & Taam, 2004), and numerical results indicate that a HG star will merge with its companion upon entering a CE (Taam & Sandquist, 2000). In the reference case, all CE with HG donors result in a merger.
We set the maximum NS mass to 2.5 , where NS that become too massive collapse as black holes. See Belczynski et al. (2008) for additional parameters dealing with the calculation of winds, tides, and stellar structure in the reference model.
2.1.2 Varied Parameters
|Model||KickaaThe kick velocity applied to a NS born in an Fe core-collapse SN is drawn randomly from a Maxwellian distribution with this dispersion velocity.||NotesbbHCE: NS accrete hypercritically in the CE; HG in CE: HG stars are allowed to survive a CE.||log()||log()||rankccThe rank of each model when ordered by .|
|6||0.5||300||HG in CE||-22.8||-23.4||16|
|7||0.5||300||HCE, HG in CE||-22.8||-23.3||11|
Note. – These models all have reference values (found in §2.1.1) for parameters unlisted here.
We explore the effects of different Fe core-collapse supernova kick velocities by testing Maxwellian distributions with smaller dispersion velocities of 150 km/s, 50 km/s, and a model with no kicks. We also test models in which the ECS kick velocity is varied, including a model in which all SN kicks (Fe core-collapse and ECS) have a velocity dispersion of 50 km/s. We test models with the following ranges for which form a partially degenerate CO core: 1.3-2.25, 1.7-1.9, 1.8-2.1, 1.83-2.25, 1.83-2.5, 2.0-2.5, 2.2-3.0, 2.5-2.7 and 1.66-3.24 . These ranges of masses are varied to include lower, higher, and wider ranges than the reference range. In dealing with the CE, we explore the parameter space of , testing both higher and lower values: 0.01, 0.05, 0.1, 0.2, 0.25, 0.3, and 1.0. We do not test values of 1.0 corresponding to an efficiency greater than 100% as this requires an additional energy source. We also test an alternative prescription, formulated by Nelemans et al. (2000) in which angular momentum, not energy, is parametrized (see also Nelemans & Tout, 2005). We further test models in which we attempt to calculate the mass accreted by a NS in a CE, set to be half the Bondi-Hoyle accretion rate. Finally, we test models allowing HG donors to survive a CE.
The maximum NS mass is set at 2.5 , but we test models lowering the limit to 2.0 , which becomes relevant in models allowing hypercritical accretion in the CE. We alter the fraction of mass lost in a stable semi-detached system by a non-compact accretor between 0, 0.1, 0.3, 1.0, and the default value of 0.5. Other variations from the reference model include decreasing the wind strengths, decreasing the effectiveness of tides, altering the mass below which the helium envelope becomes convective, using an alternative initial mass function for the secondary to create binaries with similar initial masses, altering the mass-radius relation for helium stars, increasing the specific angular momentum of matter, and various combinations of these parameters.
We also adopt a couple more assumptions which we do not vary in our set of models. One such variable is the NS birth mass (1.35 for Fe core-collapse SN and 1.26 M for an ECS). We make this assumption based on observed masses of NS, and we expect that slightly altering the birth mass of NS will not significantly affect the orbital periods and eccentricities of the DNS produced (see Özel et al., 2012, and references therein). We also do not vary the star formation rate. Although this is a simple approximation, the star formation rate probably does not vary by much more than a factor of 2 over the past 10 Gyr (Rocha-Pinto et al., 2000).
3. Model Analysis Methodology
We analyze our simulation results in the context of the DNS orbital period-eccentricity plane, and their evolutionary constraints. In Figure 1 the distribution of the entire DNS population for our reference model is shown in the first panel. The eight observed systems are overlaid as diamonds on top of the distribution. The 68.3%, 95.45%, 99.7%, and 99.994% confidence levels of the DNS population are displayed by contours in Figure 1. The characteristic band, running from low eccentricity to high eccentricity is due to the finite impulse of energy and angular momentum by a kick imparted to the second NS at its birth. Inspiral and circularization due to the emission of gravitational radiation affect systems with short orbital periods and large eccentricities. This band of systems qualitatively matches the distributions found by previous authors (Portegies Zwart & Yungelson, 1998; Belczyński & Bulik, 1999; Dewi, 2010; Kiel et al., 2010). We expand on the work of previous authors by including the three additional plots in Figure 1. These plots break down the entire population into three parts based on their evolutionary channels described in the §3.1.
Table 2 provides the model parameters for the nine models we discuss in detail hereafter. The additional 146 models are provided in Table 3. As discussed in what follows, we find that the variables that have the largest impact on our results are the He-core mass range within which a partially degenerate CO core is formed, kick velocity dispersion, whether a HG star is allowed to survive a CE, the value of , and whether hypercritical accretion is allowed in the CE.
3.1. DNS Evolutionary Channels
In the standard DNS formation model, the most important phases of evolution occur after the primary becomes a NS. Once the secondary evolves beyond the main sequence, it may enter a phase of dynamically unstable mass transfer. As mentioned in §1, we completely ignore non-interacting systems because the short orbital periods in observed DNS require mass transfer. The secondary will evolve onto the helium main sequence, possibly filling its Roche lobe upon evolving off the helium main sequence in the so-called Case BB mass transfer (Delgado & Thomas, 1981). After forming a degenerate core, the star will collapse forming a NS in either an Fe core-collapse SN or an ECS.
Variations from the standard paradigm have been discussed in the literature or seen in population synthesis studies (Bhattacharya & van den Heuvel, 1991; Brown, 1995; Portegies Zwart & Yungelson, 1998; Belczynski et al., 2002). We find the scenario outlined by Brown (1995) in which the initial mass ratio is close to unity, and the individual stars in the binary evolve off the main sequence together, forming double common envelopes occurs in less than 1 percent of systems, and ignore them in our analysis, (however see Schwab et al., 2010). While we find some NS form through the accretion induced collapse of a WD the population is insignificant and we disregard them. Taking the standard formation paradigm outlined above, we can break down our models into three different evolutionary tracks depending on whether or not the system went through stable mass transfer after a CE (Case BB mass transfer), and, if it did, whether the second NS was formed in an ECS or an Fe core-collapse SN.
Channel I: As the simplest of our three formation channels, the systems that are created this way have the fewest constraints on their evolution. Therefore this DNS population covers the widest range in orbital period and eccentricity which can be seen in the second panel in Figure 1. Following the formation of the first NS, the secondary evolves off the main sequence; the orbits of this population are wide enough that the massive star fills its Roche lobe when a significant convective envelope has developed and a CE ensues. The NS companion is massive enough that it eventually forms an Fe core before collapsing into the second NS in the system.
Channel II: In this channel, the binary orbit is tight enough following the CE event that once core He burning has finished the star expands, filling its Roche lobe again in Case BB mass transfer; provided the He star has not developed a fully convective envelope, mass transfer onto the first NS is stable (see Ivanova et al., 2003). The secondary evolves to form an Fe core and eventually collapses in a SN. Comparing panel two to three in Figure 1 shows that this stable mass transfer phase produces a DNS population with shorter orbital periods than that of Channel I.
Channel III: These systems follow a similar evolution to Channel II, except prior to core-collapse, the secondaries have lower He core masses, leading to NS formation through an ECS. These lower mass He star progenitors always expand significantly leading to the second phase of mass transfer onto the primary NS. Unless the system has a very small orbital separation, this phase of mass transfer is stable (Dewi et al., 2002; Ivanova et al., 2003). Since ECS kick velocities are small, the resulting DNS orbital geometries are more dependent on the mass loss of the secondary during NS formation. Due to the tight constraints on the mass of the secondary immediately prior to the SN (the secondary star has a mass within a narrow range such that carbon is allowed to burn, but oxygen is not), the resulting population has a very narrow distribution of geometries. This can be seen by the relatively narrow distribution of systems in Figure 1. Combined with the small kick velocities, systems going through this evolutionary channel tend to have systematically lower eccentricities.
3.2. Evolutionary Constraints
A few known DNS systems have such characteristics that a number of different groups have been able to derive firm constraints on their evolutionary history, and we are able to confidently place them in the channels described here. Having this knowledge allows us to take a step further and, for the first time, impose evolutionary constraints on the population formation. In what follows we summarize the origin of these evolutionary constraints on the three best studied DNS systems. Once these constraints are imposed, we use the best models to draw conclusions about the most likely channels through which the rest of the DNS systems have evolved (§4.3).
J07373039: The binary pulsar is probably the most widely studied DNS in terms of its evolutionary history, which is now generally agreed upon. After the primary star became a NS, the secondary had to evolve off its main sequence filling its Roche lobe in the process and forming a CE. Both Dewi & van den Heuvel (2004) and Willems et al. (2004) conclude, based on the current geometry of the system, that it must have filled its Roche lobe prior to the second SN. Based on the small proper motion measurement and close proximity to the galactic disk, Piran & Shaviv (2005) conclude that the He star progenitor of J07373039B was probably , and the system received a SN kick km/s at birth. This possibility, although more complex in its origin than claimed by Piran & Shaviv (2005), was confirmed by other studies (Willems et al., 2006; Stairs et al., 2006; van den Heuvel, 2007; Wong et al., 2010; Ferdman et al., 2013; Dall’Osso et al., 2014). Limits on the He star progenitor mass for this system preclude the second born NS from having ever formed an Fe core (Nomoto, 1987); prior to core-collapse, it was only massive enough to form an ONe core. This constraint combined with its very low natal kick requires that this system formed in an ECS. Recently, a population synthesis study by Dewi (2010) showed that the He star progenitor of J07373039B could not have overflowed its Roche lobe in a second CE event; instead the mass transfer phase by the He progenitor must have been stable. Thus, we can conclude that J07373039 was formed through Channel III.
B1534+12: This system, too, has been studied quite extensively (Willems et al., 2004; Thorsett et al., 2005; Stairs et al., 2006; Wong et al., 2010). All groups agree that the only possible way to explain all the system parameters was if the secondary had gone through a phase of stable mass transfer as a He star. However, this system differs from J0737+3039 in two ways. First, the second NS has a mass of 1.35 , consistent with the expected remnant of an Fe core-collapse SN. Second, the natal kick magnitude for this second NS is constrained to be high, on order of several hundred km/s. Based on our assumptions about ECS, we place B1534+12 within Channel II.
B1913+16: Although many authors have commented on the restrictions of the kick velocity and direction required to form B1913+16, the evolutionary constraints on the system’s formation channel are less strong (Fryer & Kalogera, 1997; Willems et al., 2004; Ihm et al., 2006; Wong et al., 2010). According to these studies, the kick velocity forming the second NS must have been large ( km/s), and the mass of the secondary NS is the expected result of an Fe core-collapse SN, so an ECS is ruled out for the formation of the secondary NS. However, since the mass transfer history of this system has not yet been constrained, we can only restrict B1913+16 to have evolved through either Channel I or Channel II.
3.3. Statistical Analysis
To quantitatively analyze the comparative goodness of fit for each model, we adapt the Bayesian analysis used in Ihm et al. (2006) to fit for two independent parameters, orbital period and eccentricity. We begin by applying Bayes’s Theorem to our results:
where is the probability of the model being correct given the data, is the probability of the data given the model, is the prior probability of the model, and is a normalizing constant that is independent of model, . Our analysis, however, includes prior evolutionary constraints, so Bayes’s Theorem becomes:
where is the value which we are calculating: the probability that our model is correct, given the data. is the probability of the observed values given our model. is the probability that, given our model, the evolutionary constraints for the data are satisfied. is our prior, the a priori probability that each particular model is correct. We give each model the same prior probability. Let
Because we have given each model the same prior probability, is independent of model. Let be
Since is independent of model, the values of give the relative probabilities of each model. Now, we substitute orbital period, , and eccentricity, , for our data in our equation for . Due to the independence of each observed DNS, the probability of the data set given the model is equal to the product of the probabilities of each independent system:
It is necessary to include a subscript with because there are different evolutionary constraints on each observed system. can be easily determined as the fraction of systems that go into each evolutionary channel.
The errors for the values of eccentricity and orbital period are extremely small and can be ignored. Therefore, calculating reduces to finding the probability density of our model distributions at each point in the two dimensional space corresponding to the observed DNS. Because our data space is two-dimensional, and each model carries no additional parameters, we do not need to use MCMC-type methods in the computation of the model posterior; instead we calculate by binning the simulated systems into two dimensions for each evolutionary channel. We use a bin size of in eccentricity and in orbital period. Varying the bin size gives a good handle on the inherent error of our statistical method. Using equation (5), the value we are calculating, , is equal to the product of the probability density of each bin that corresponds to the observed systems. We find that, although our models have values of that vary by several orders of magnitude, changes within a single model by at most a factor of a few as the bin size varies over reasonable ranges, so our model rankings are insensitive to our choice of binning.
It is useful to determine the dependencies of our results on our use of evolutionary constraints. To do this, we marginalize over the nuisance parameter of the evolutionary channels. The analysis is performed similarly, the only exception being our entire model distribution is binned in two dimensions without evolutionary constraints. In this case, can be defined as:
The values for (hereafter ) will, in general, differ from the values we find for (hereafter ).
In what follows we analyze our results and draw conclusions in the context of a number of different questions. First we examine how the model behavior on the plane is affected by a few critical model parameters. We then use the likelihood calculation for two cases, one that includes only the measurements as constraints and one that accounts for the evolutionary constraints of the three systems discussed in § 3.2, and we identify the most favored models and discuss their parameters. For these top models, we also calculate the branching ratios through the three evolutionary channels and try to identify the most likely channel for all eight DNS in the known sample.
4.1. Influence of Key Model Parameters
The models in Figures 2 and 3 demonstrate the influence of the energy efficiency in CE evolution, SN kicks, as well as hypercritical NS accretion and survival through the CE phase when the donor star is in the HG.
As expected, tighter orbits are produced by models with less efficient CE phases, those with smaller values. This is reflected in Models 2, 3, and 4 with , , and , respectively in Figure 2. This figure further demonstrates the effect on branching ratios between the models; smaller CE efficiencies lead to proportionately fewer systems becoming DNS through Channel III. The interplay of these two effects demonstrate the difficulty in finding models that reproduce J07373039. Models 4 and 5, which also have difficulty of producing systems similar to J07373039, nonetheless remain viable models due to their increased ability to reproduce other DNSs, particularly J15184904 and J18292456.
Smaller NS natal kick velocities result preferably in less eccentric DNS orbits. Model 8 with a kick velocity dispersion of 50 km/s in Figure 3 shows how the highest concentrations of DNS tend to have eccentricities smaller than 0.5. Although these low kick models better reproduce the small eccentricity DNS, they have difficulty reproducing B1913+16. Low kick models alter the branching ratios, creating more DNS through Channel I. Model 9 in Figure 3, which applies kicks with a 150 km/s dispersion, shows an intermediate between Model 8 and the standard model.
Model 5 in Figure 2 demonstrates the shift toward slightly longer orbital periods when NS are allowed to accrete hypercritically in a CE. Such accretion increases the NS mass, reduces the envelope mass to be ejected and therefore requires a smaller degree of orbital contraction for envelope ejection. On the other hand, allowing HG stars to survive the CE phase, if there is enough orbital energy initially, has an almost negligible effect as seen by Model 6 in Figure 3. There is a slightly higher density of systems at the shortest orbital periods in Channel III. Model 7, which both allows NS to accrete hypercritically as well as allows HG stars survive a CE, shows the combination of these two effects, as expected.
In their study of double compact object merger rates, Dominik et al. (2012) find that only relatively close binaries will overfill their Roche lobe while the donor is on the HG. If these donors are allowed to survive a CE, they will produce tight binaries that merge relatively quickly due to gravitational wave radiation (as evidenced by their short calculated delay time distribution), resulting in an increased DNS merger rate. Since these systems merge soon after their formation, they are unlikely to be observed as DNS. Therefore, while this evolutionary channel may significantly contribute to the DNS merger rate, its effect on the distribution of observed DNS in space is negligible.
Using a range for of 1.83-2.5 populates Channel III, while the number of systems formed through Channel III drops precipitously for a slightly lower range of 1.83-2.25 . We further find that while we test several different ranges, we find that only those models with a range for of 1.83-2.5 and 2.0-2.5 produce DNS through Channel III in reasonable numbers. If the He star mass range is not high enough, either the partially degenerate CO core will not become massive enough (1.08 ) to begin stable C burning, or the resulting ONe core will not reach the Chandrasekhar limit (1.38 ) and leave an ONe WD remnant.
Constraints on were independently proposed by Linden et al. (2009) who analyzed the population of high mass X-ray binaries in the Small Magellanic Cloud. They found that if X-ray binaries were not formed through an ECS, the numbers of observed systems with ages 20-60 Myr could not be explained. However, when using a range for that allowed X-ray binaries to be formed through ECS, X-ray binaries were formed with the correct age distribution. Although they did not determine the best range for , their fiducial range of 1.83-2.25 compares to our preferred ranges for of 1.83-2.5 and 2.0-2.5 .
4.2. Best-fit Models and their Parameters
Following the analysis of §3.3, we evaluate each of our 155 models by calculating the two likelihoods, and , without and with the evolutionary constraints, respectively; the values are shown in Figure 4 in decreasing order across model number. The ordering is performed separately in each panel. When evolutionary constraints are included in the analysis, the number of viable models drops from 100 to 70. We select these 70 models as their value is within three orders of magnitude of the highest likelihood model. Our remaining discussion will be confined to these “best-fit” models indicated by the dotted line in the bottom panel of Figure 4. Table 2 shows the input parameters of the top 70 models, including the reference model.
Within the top 70 models, there is little preference either for or against models allowing hypercritical accretion in the CE or allowing HG stars to survive a CE. Furthermore, the secondary parameters have little effect on the resulting distributions. These include changing the mass beyond which a helium core is fully convective, using an angular momentum prescription to calculate the CE evolution, altering the efficiency of mass accretion during stable mass transfer, altering the wind mass loss prescription and decreasing the maximum allowed neutron star mass.
Two different ranges for dominate the top models in Table 2: 1.83-2.5 and 2.0-2.5 . This is due to the difficulty for other ranges to create DNS through Channel III, a necessity for forming J07373039.
There is an interplay between and the natal kick velocities applied in our models: low kick models push the distributions toward smaller eccentricities, while high CE efficiency models push the distributions toward longer orbital periods. Combined, these models make it difficult to form J07373039 through Channel III. Therefore, when evolutionary constraints are applied to the models, the models most often eliminated are those with and a low kick velocity. This is consistent with independent but similar constraints derived by Linden et al. (2009) who argue that large Fe core-collapse kicks and small ECS kicks are required to explain the population of high mass X-ray binaries within the Small Magellanic Cloud. Nevertheless, there are viable models with all three kick velocities that we test. We note that in models with a kick dispersion velocity of 300 km/s, the bulk of DNS are formed with , while only two of the eight known DNS have such high eccentricities. We show in §5.2 that this is unlikely to be due to a bias against detecting high eccentricity DNS.
Figure 5 shows the branching ratios of the top models, ordered by decreasing . The majority of DNS in these top models are formed through channel II shown as a green line. These systems are characterized by partially recycled primary NS, and secondary NS with masses of 1.35 , consistent with many of the DNS in Table 1. The red line in Figure 5 shows that in these models, DNS are formed in appreciable numbers through Channel III.
4.3. Most Likely Evolutionary Channels for DNS Systems
One motivation for the present study is to examine whether we can identify the statistically-favored evolutionary models for the DNS systems, for which, at present, it is impossible to identify their evolutionary history based on their individual studies. With this in mind, we re-assess the branching ratios for the best-fit models (which have been “calibrated” against the systems with known evolutionary channels) in direct connection with the observed systems. For each observed system and for each best-fit model, we select the sub-population of model systems with eccentricities and orbital periods in the same 2D bin as is used to calculate our values in §3.3. The evolutionary channel branching ratios are calculated as the proportion of systems falling within this bin that go through each evolutionary channel. Varying the bin size affects the evolutionary channel branching ratios at the 10-15% level.
The results for all 70 viable models are shown in Figure 6. For systems such as B1913+16 and J1518+4904, one channel is clearly favored statistically, although this is not necessary the channel that actually formed each system. A few of these models are peculiar. The models at 31 and 48 are models that calculate CE evolution based on the conservation of angular momentum not energy. Although our statistical analysis gives each model an equal prior, modern discussions generally disfavor this prescription (Webbink, 2008; Woods et al., 2011; Ivanova et al., 2013).
J07373039: This is the prime example illustrating that the statistically favored channel is not necessarily the true channel. In this case, the former is Channel II, while the latter is Channel III, based on the analysis of all its observed characteristics (as discussed in § 3.2). In Figure 6, the red line representing Channel III shows that in all models, no more than half of all DNS are formed through this pathway.
B1534+12: The requirement of B1534+12 to form through Channel II is in agreement with the branching ratio results shown in Figure 6; in this case the statistically favored channel is also the true channel. As discussed in § 3.2 the mass measurements and NS kick constraints indicate that B1534+12 could not have been formed in an ECS, precluding formation through Channel III, the other likely evolutionary channel in Figure 6.
B1913+16: Prior evolutionary constraints restrict B1913+16 to be created through either Channel I or II. In Figure 6 we see that, of these two, Channel II is greatly favored from a statistical point of view. The models at 31 and 48 correspond to models that calculate CE evolution based on an angular momentum conservation description. Hence we can say that B1913+16 most likely went through a phase of stable mass transfer prior to the second SN event. This can be attributed to both Channel II having a much higher branching ratio as well as the systems going through Channel I being formed at larger orbital periods.
J17562251: Recent observations by Ferdman et al. (2014) indicate the pulsar companion in J17562251 has a mass of 1.23 . The low mass combined with their newly measured tangential velocity of 20 km s suggests that the pulsar companion may have been born in an ECS. If that is the case, J17562251 must have been formed through Channel III. The branching ratio for Channel III (red in Figure 6) varies between 10% and 60% for J17562251, indicating this possibility. If we add the constraint requiring J17562251 to be formed through Channel III, the resulting values do not change significantly.
J1906+0746: J1906+0746 is unique in that the observed pulsar is the second-born NS in the system. A mass measurement for the pulsar indicates it may have been born in an ECS, constraining the system to Channel III. Figure 6 indicates that this system almost certainly went through either Channel II or III. Through future analysis of the galactic position and velocity of J1906+0746, it may be possible to determine if it went through a phase of stable mass transfer, and therefore Channels II or III. However, with the available information we cannot differentiate between the two channels, as the branching ratios in Figure 6 indicate both are possible.
J18111736: The branching ratios in Figure 6 indicate that J18111736 likely evolved through Channel I or II. Its relatively large orbital period means formation through Channel III is rare for most models. Models 7, 8, and 9 in Figure 3 show this is because of the difficulty of forming a DNS with large orbital periods and eccentricities through the low velocity kick applied in an ECS. The models that allow for the formation of J18111736 through Channel III tend to be those with higher kick velocities, and hence larger ECS kicks. Based on their measurements of the relativistic periastron advance and the low derived system mass, Corongiu et al. (2007) suggest that J18111736 was born with a low velocity kick. Although uncertainties in the current mass measurements for the system are too large to determine evolutionary constraints, Corongiu et al. (2007) calculate that a second post-Keplerian parameter may be measureable in the near future. If future observations indicate this system was formed through Channel III, it could provide an important constraint on ECS kicks.
J1518+4904: The 8.6 day orbital period and eccentricity of 0.249 for J1518+4904 make its formation difficult. The distributions in Figures 2 and 3 show that most DNS are formed at either smaller orbital periods or larger eccentricities. Although a few models form J1518+4904 through Channel II or III, Figure 6 indicates it was most likely formed through Channel I. Observational errors on the mass estimates are too large to indicate whether the companion was born in an ECS event or not.
J1829+2456: Depending on the individual model, Figure 6 shows that J1829+2456 formed through any of the three evolutionary channels. Poor mass constraints do not allow any further indication. Future observations may provide further constraints on its evolution.
5.1. Partial Recycling
While there is little observational evidence to distinguish between different evolutionary channels, one such piece of evidence could be the degree to which the first-born NS has been recycled. Ivanova (2011) has argued that the mass losing star in a CE will go through a phase of thermal readjustment upon envelope expulsion. During this phase a NS companion can be mildly recycled. Investigation of Table 1 shows that there may be in fact a bifurcation spin periods in the primary NS: those with 50 ms and those with 100 ms. If this argument is accurate, then J18111736 went through Channel I, avoiding any stable mass transfer after the CE phase, while the other systems with shorter spin periods went through Channel II or III. J1906+0746 is an exception because we do not know the spin period of the primary NS.
Comparison with Figure 6 shows that this idea is consistent with most of the results of our analysis here. In roughly half of our viable models, Channel I is preferred for J18111736. With a few exceptions, other systems prefer Channel II (in agreement with their shorter spin periods) except for J1518+4904 which our results here indicate was likely formed through Channel I. Either this system, too, formed through Channel II, or the primary NS in J1518+4904 accreted enough material in a CE to recycle the primary to its current spin period of 40.9 ms.
5.2. Observational Biases
Our statistical analysis discussed in §3.3 relies on the assumption that all DNS have an equal detection likelihood. For example, a strong bias against detection of DNS with large eccentricities means that high kick velocity models, which produce on average more eccentric DNS, will not be as disfavored as the values indicate. In principle our DNS distributions could be convolved with any observational bias for a more accurate comparison to the observed DNS.
One such bias deals with the observability of pulsars in DNS which is a function of a number of intrinsic parameters such as the beaming fraction, magnetic field strength, age, and spin period as well as particular parameters such as the distance and direction from the Sun. However, since we ignore non-interacting DNS, our model assumes that every DNS contains a partially recycled pulsar, each with an equal probability of being observed. Absent any deeper understanding of the connection between a pulsar’s characteristics and its prior evolution, we take this to be a reasonable approximation.
An additional bias could be caused by Doppler smearing of pulsars in short orbital periods. Due to the high accelerations in such systems over the course of a single observation, close binaries suffer a decreased detection efficiency. Bagchi et al. (2013) analyzed the detectability of binary pulsars in surveys, finding that modern acceleration searches (e.g. PRESTO) were greater than 80% efficient at detecting pulsars in DNS for typical parameters. Interestingly, they found systems with the shortest orbital periods and lowest eccentricities suffered the strongest biases, which still have a detection efficiency greater than 50%.
A potentially more pernicious bias could be introduced by the limited sky coverage of pulsar surveys, typically close to the Galactic plane. Binaries that gain large systemic velocities due to the NS natal kicks may therefore be missed. Such DNS would likely have significant eccentricity and a large orbital period, possibly helping to explain the relative dearth of such systems despite the predictions by our simulations. We first investigate the systemic velocities of our resulting DNS populations.
The distribution of systemic velocities, normalized to the highest bin, for Models 8 ( km/s), 9 ( km/s), and 1 ( km/s) are shown in Figure 7. We expect these three models to be representative of the different kick velocity models of our entire set. Calculated from the equations in Kalogera (1996), this systemic velocity arises from both the mass lost from the collapsing object as well as the kick applied to the newly born NS. Model 8 produces DNS with systemic velocities peaking at 35 km/s. Model 9 and 1 both produce DNS with systemic velocities peaking at 75 km/s, although the higher kick velocity Model 1 has a high velocity tail. Can these systemic velocities create a bias against detection strong enough to significantly affect our quantitative results?
To determine this magnitude of any effect, we create a model Milky Way-like galaxy composed of a three-component potential defined in Appendix A of Wex et al. (2000). DNS are born in a double exponential, axisymmetric disk with a scale length of 2.80 kpc and scale height of 70 pc, and are given an initial velocity corresponding to circular rotation. The motion of the DNS through the Galaxy is calculated using a Runge-Kutta fourth-order scheme. The binary is evolved for 20 Myr (the approximate evolution time of the binary prior to the birth of the second NS) before a systemic velocity is applied. The systemic velocity is applied to the system in a randomly chosen direction. The system is then evolved for a Gyr and the location of the binary recorded throughout this evolution at intervals of a Myr. This location is translated into a galactic longitude and latitude as referenced from the Sun, placed at 8.5 kpc from the center of the Galaxy. We repeat this process, simulating the motion of 1000 binaries throughout the Galaxy, each with a different randomly-chosen initial position and systemic velocity direction. We then determine the percentage of the resultant Galactic longitudes and latitudes that fall within a region near the Galactic Plane, using as a fiducial region the Parkes multibeam pulsar survey area (Manchester et al., 2001): and . Figure 8 shows this percentage over the relevant range of systemic velocities. We find that when no velocity is applied, a DNS has a 93% chance of falling within the Parkes multibeam pulsar survey region. That number falls to 85% for systemic velocities of 100 km/s and 69% for systemic velocities of 200 km/s.
We test the strength of this bias on Models 8, 9, and 1 with Maxwellian kick dispersion velocities of 50, 150, and 300 km/s, respectively. For each simulated system in the three models, we find the detection probability from Figure 8 corresponding to its individual systemic velocity. We then convolve observational biases into our distributions in space by selecting a subset of our original distribution. Each simulated DNS has a chance equal to its calculated detection probability of being included into the subset. In Model 1 (with the highest kick velocities and hence the strongest expected effects of the three models) roughly 80% of all systems are unaffected. The effects of this observational bias are seen in Figure 9; the top panels provide the full distributions of DNS for our three test models while the bottom panels show the distributions convoluted with the observational bias. The resulting distributions are nearly qualitatively identical. Since our quantitative results are robust to a factor of several orders of magnitude, we can safely ignore this observational bias as potentially altering our conclusions in this study.
We generate populations of DNS using a binary population synthesis code. In total, we test 155 models. The simulated distributions are broadly consistent with the orbital periods and eccentricities of the known systems. These distributions are split into three separate evolutionary channels based on two criteria: if the secondary star entered a phase of stable Roche lobe overflow as a He star and if the secondary formed a NS in an Fe core-collapse SN or an ECS. Previous authors have constrained the evolutionary histories of three of the eight known DNS: B1913+16, B1534+12, and J07373039.
For the first time, we combine evolutionary constraints on the three known systems with a Bayesian analysis to quantitatively estimate the relative likelihood of each model. Of our 155 total models, we find that 70 can create B1913+16, B1534+12, and J07373039 through their constrained evolutionary channels.
We use these 70 viable models to further constrain binary evolution parameters. We find that a common envelope efficiency is effectively ruled out. Although values between 0.3 and 1.0 are allowed, 0.5 tends to fit the data best. Since evolutionary codes typically find for the high-mass AGB donor stars within a CE, DNS formation may require efficient CE evolution. More work is required to determine if this requires greater than unity. The most constraining parameter is , the He core mass which allows the formation of a partially degenerate CO core. Of all the mass ranges we test, only 1.83-2.5 and 2.0-2.5 allow for the formation of J07373039 in reasonable numbers. We test models with NS natal kicks drawn from a Maxwellian distribution. Models survive with all three tested dispersion velocities (50, 150, and 300 km/s).
We can further use our 70 viable models to constrain which evolutionary channels formed each of the eight known DNS. We find that B1534+12 likely went through Channel II, in agreement with its evolutionary constraint. Although previous work constrains B1913+16 to have gone through either Channel I or II, our results here indicate it most likely went through Case BB mass transfer in our Channel II. Although difficult to form in general, J0737-3039 could have been formed through Channel III, in agreement with its evolutionary constraint. J17562251 and J1906+0746 may have been formed through Channel III, consistent with the idea that ECS produce low mass NS. J18111736 likely went through Channel I or II, however if future observations and analysis indicate it was born through Channel III, it could place a strong constraint on ECS kick velocities. Our models indicate J1518+4904 most likely formed through Channel I. It is clear that understanding the evolutionary history of known DNS systems provides strong constraints on population models and helps us further constrain NS formation.
As already discussed J07373039 appears to have formed through channel III, which systematically has a lower branching ratio compared to channel II (no more than about 50% as shown in Figure 6). These low branching ratio values may seem in contrast with early indications that the formation rates of pulsar binaries similar to J07373039 dominate the total DNS by a factor of 6-7 (Kalogera et al., 2004). However, the most recent analysis of the double pulsar system (Kim et al., 2013) concludes that current beaming constraints lead to the empirical DNS rates being comparable for J07373039 and B1913+16.
- Bagchi et al. (2013) Bagchi, M., Lorimer, D. R., & Wolfe, S. 2013, MNRAS, 432, 1303
- Belczyński & Bulik (1999) Belczyński, K. & Bulik, T. 1999, A&A, 346, 91
- Belczynski et al. (2010) Belczynski, K., Bulik, T., Fryer, C. L., Ruiter, A., Valsecchi, F., Vink, J. S., & Hurley, J. R. 2010, ApJ, 714, 1217
- Belczynski et al. (2002) Belczynski, K., Kalogera, V., & Bulik, T. 2002, ApJ, 572, 407
- Belczynski et al. (2008) Belczynski, K., Kalogera, V., Rasio, F. A., Taam, R. E., Zezas, A., Bulik, T., Maccarone, T. J., & Ivanova, N. 2008, ApJS, 174, 223
- Bhattacharya & van den Heuvel (1991) Bhattacharya, D. & van den Heuvel, E. P. J. 1991, Phys. Rep., 203, 1
- Brown (1995) Brown, G. E. 1995, ApJ, 440, 270
- Burgay et al. (2003) Burgay, M., D’Amico, N., Possenti, A., Manchester, R. N., Lyne, A. G., Joshi, B. C., McLaughlin, M. A., Kramer, M., Sarkissian, J. M., Camilo, F., Kalogera, V., Kim, C., & Lorimer, D. R. 2003, Nature, 426, 531
- Champion et al. (2004) Champion, D. J., Lorimer, D. R., McLaughlin, M. A., Cordes, J. M., Arzoumanian, Z., Weisberg, J. M., & Taylor, J. H. 2004, MNRAS, 350, L61
- Corongiu et al. (2004) Corongiu, A., Kramer, M., Lyne, A. G., Löhmer, O., D’Amico, N., & Possenti, A. 2004, Memorie della Societa Astronomica Italiana Supplementi, 5, 188
- Corongiu et al. (2007) Corongiu, A., Kramer, M., Stappers, B. W., Lyne, A. G., Jessner, A., Possenti, A., D’Amico, N., & Löhmer, O. 2007, A&A, 462, 703
- Dall’Osso et al. (2014) Dall’Osso, S., Piran, T., & Shaviv, N. 2014, MNRAS, 438, 1005
- de Kool (1990) de Kool, M. 1990, ApJ, 358, 189
- De Marco et al. (2011) De Marco, O., Passy, J.-C., Moe, M., Herwig, F., Mac Low, M.-M., & Paxton, B. 2011, MNRAS, 411, 2277
- Delgado & Thomas (1981) Delgado, A. J. & Thomas, H.-C. 1981, A&A, 96, 142
- Deloye & Taam (2010) Deloye, C. J. & Taam, R. E. 2010, ApJ, 719, L28
- Dewi (2010) Dewi, J. D. M. 2010, New A, 54, 145
- Dewi et al. (2002) Dewi, J. D. M., Pols, O. R., Savonije, G. J., & van den Heuvel, E. P. J. 2002, MNRAS, 331, 1027
- Dewi & van den Heuvel (2004) Dewi, J. D. M. & van den Heuvel, E. P. J. 2004, MNRAS, 349, 169
- Dominik et al. (2012) Dominik, M., Belczynski, K., Fryer, C., Holz, D. E., Berti, E., Bulik, T., Mandel, I., & O’Shaughnessy, R. 2012, ApJ, 759, 52
- Ferdman et al. (2013) Ferdman, R. D., Stairs, I. H., Kramer, M., Breton, R. P., McLaughlin, M. A., Freire, P. C. C., Possenti, A., Stappers, B. W., Kaspi, V. M., Manchester, R. N., & Lyne, A. G. 2013, ArXiv e-prints
- Ferdman et al. (2014) Ferdman, R. D., Stairs, I. H., Kramer, M., Janssen, G. H., Bassa, C. G., Stappers, B. W., Demorest, P. B., Cognard, I., Desvignes, G., Theureau, G., Burgay, M., Lyne, A. G., Manchester, R. N., & Possenti, A. 2014, MNRAS, 443, 2183
- Flannery & van den Heuvel (1975) Flannery, B. P. & van den Heuvel, E. P. J. 1975, A&A, 39, 61
- Fryer & Kalogera (1997) Fryer, C. & Kalogera, V. 1997, ApJ, 489, 244
- Fryer et al. (2012) Fryer, C. L., Belczynski, K., Wiktorowicz, G., Dominik, M., Kalogera, V., & Holz, D. E. 2012, ApJ, 749, 91
- Ge et al. (2010) Ge, H., Hjellming, M. S., Webbink, R. F., Chen, X., & Han, Z. 2010, ApJ, 717, 724
- Hobbs et al. (2005) Hobbs, G., Lorimer, D. R., Lyne, A. G., & Kramer, M. 2005, MNRAS, 360, 974
- Hulse & Taylor (1975) Hulse, R. A. & Taylor, J. H. 1975, ApJ, 195, L51
- Hurley et al. (2000) Hurley, J. R., Pols, O. R., & Tout, C. A. 2000, MNRAS, 315, 543
- Iben & Livio (1993) Iben, J. I. & Livio, M. 1993, PASP, 105, 1373
- Ihm et al. (2006) Ihm, C. M., Kalogera, V., & Belczynski, K. 2006, ApJ, 652, 540
- Ivanova (2011) Ivanova, N. 2011, ApJ, 730, 76
- Ivanova et al. (2003) Ivanova, N., Belczynski, K., Kalogera, V., Rasio, F. A., & Taam, R. E. 2003, ApJ, 592, 475
- Ivanova & Chaichenets (2011) Ivanova, N. & Chaichenets, S. 2011, ApJ, 731, L36
- Ivanova et al. (2013) Ivanova, N., Justham, S., Chen, X., De Marco, O., Fryer, C. L., Gaburov, E., Ge, H., Glebbeek, E., Han, Z., Li, X.-D., Lu, G., Marsh, T., Podsiadlowski, P., Potter, A., Soker, N., Taam, R., Tauris, T. M., van den Heuvel, E. P. J., & Webbink, R. F. 2013, A&A Rev., 21, 59
- Ivanova & Taam (2004) Ivanova, N. & Taam, R. E. 2004, ApJ, 601, 1058
- Jacoby et al. (2006) Jacoby, B. A., Cameron, P. B., Jenet, F. A., Anderson, S. B., Murty, R. N., & Kulkarni, S. R. 2006, ApJ, 644, L113
- Janssen et al. (2008) Janssen, G. H., Stappers, B. W., Kramer, M., Nice, D. J., Jessner, A., Cognard, I., & Purver, M. B. 2008, A&A, 490, 753
- Kalogera (1996) Kalogera, V. 1996, ApJ, 471, 352
- Kalogera (2000) —. 2000, ApJ, 541, 319
- Kalogera et al. (2004) Kalogera, V., Kim, C., Lorimer, D. R., Burgay, M., D’Amico, N., Possenti, A., Manchester, R. N., Lyne, A. G., Joshi, B. C., McLaughlin, M. A., Kramer, M., Sarkissian, J. M., & Camilo, F. 2004, ApJ, 601, L179
- Kasian (2008) Kasian, L. 2008, in American Institute of Physics Conference Series, Vol. 983, 40 Years of Pulsars: Millisecond Pulsars, Magnetars and More, ed. C. Bassa, Z. Wang, A. Cumming, & V. M. Kaspi, 485–487
- Keith et al. (2009) Keith, M. J., Kramer, M., Lyne, A. G., Eatough, R. P., Stairs, I. H., Possenti, A., Camilo, F., & Manchester, R. N. 2009, MNRAS, 393, 623
- Kiel et al. (2010) Kiel, P. D., Hurley, J. R., & Bailes, M. 2010, MNRAS, 406, 656
- Kiel et al. (2008) Kiel, P. D., Hurley, J. R., Bailes, M., & Murray, J. R. 2008, MNRAS, 388, 393
- Kim et al. (2013) Kim, C., Bhakthi Pranama Perera, B., & McLaughlin, M. A. 2013, ArXiv e-prints
- Kitaura et al. (2006) Kitaura, F. S., Janka, H., & Hillebrandt, W. 2006, A&A, 450, 345
- Linden et al. (2009) Linden, T., Sepinsky, J. F., Kalogera, V., & Belczynski, K. 2009, ApJ, 699, 1573
- Livio & Soker (1988) Livio, M. & Soker, N. 1988, ApJ, 329, 764
- Lorimer et al. (2006) Lorimer, D. R., Stairs, I. H., Freire, P. C., Cordes, J. M., Camilo, F., Faulkner, A. J., Lyne, A. G., Nice, D. J., Ransom, S. M., Arzoumanian, Z., Manchester, R. N., Champion, D. J., van Leeuwen, J., Mclaughlin, M. A., Ramachandran, R., Hessels, J. W., Vlemmings, W., Deshpande, A. A., Bhat, N. D., Chatterjee, S., Han, J. L., Gaensler, B. M., Kasian, L., Deneva, J. S., Reid, B., Lazio, T. J., Kaspi, V. M., Crawford, F., Lommen, A. N., Backer, D. C., Kramer, M., Stappers, B. W., Hobbs, G. B., Possenti, A., D’Amico, N., & Burgay, M. 2006, ApJ, 640, 428
- Manchester et al. (2001) Manchester, R. N., Lyne, A. G., Camilo, F., Bell, J. F., Kaspi, V. M., D’Amico, N., McKay, N. P. F., Crawford, F., Stairs, I. H., Possenti, A., Kramer, M., & Sheppard, D. C. 2001, MNRAS, 328, 17
- Miyaji & Nomoto (1987) Miyaji, S. & Nomoto, K. 1987, ApJ, 318, 307
- Miyaji et al. (1980) Miyaji, S., Nomoto, K., Yokoi, K., & Sugimoto, D. 1980, PASJ, 32, 303
- Nelemans & Tout (2005) Nelemans, G. & Tout, C. A. 2005, MNRAS, 356, 753
- Nelemans et al. (2000) Nelemans, G., Verbunt, F., Yungelson, L. R., & Portegies Zwart, S. F. 2000, A&A, 360, 1011
- Nomoto (1984) Nomoto, K. 1984, ApJ, 277, 791
- Nomoto (1987) —. 1987, ApJ, 322, 206
- Osłowski et al. (2011) Osłowski, S., Bulik, T., Gondek-Rosińska, D., & Belczyński, K. 2011, MNRAS, 413, 461
- Özel et al. (2012) Özel, F., Psaltis, D., Narayan, R., & Santos Villarreal, A. 2012, ApJ, 757, 55
- Peters (1964) Peters, P. C. 1964, Phys. Rev., 136, B1224
- Pfahl et al. (2002a) Pfahl, E., Rappaport, S., & Podsiadlowski, P. 2002a, ApJ, 573, 283
- Pfahl et al. (2002b) Pfahl, E., Rappaport, S., Podsiadlowski, P., & Spruit, H. 2002b, ApJ, 574, 364
- Piran & Shaviv (2005) Piran, T. & Shaviv, N. J. 2005, Phys. Rev. Lett., 94, 051102
- Podsiadlowski et al. (2004) Podsiadlowski, P., Langer, N., Poelarends, A. J. T., Rappaport, S., Heger, A., & Pfahl, E. 2004, ApJ, 612, 1044
- Portegies Zwart & Yungelson (1998) Portegies Zwart, S. F. & Yungelson, L. R. 1998, A&A, 332, 173
- Rocha-Pinto et al. (2000) Rocha-Pinto, H. J., Scalo, J., Maciel, W. J., & Flynn, C. 2000, A&A, 358, 869
- Scalo (1986) Scalo, J. M. 1986, Fund. Cosmic Phys., 11, 1
- Schwab et al. (2010) Schwab, J., Podsiadlowski, P., & Rappaport, S. 2010, ApJ, 719, 722
- Stairs et al. (2006) Stairs, I. H., Thorsett, S. E., Dewey, R. J., Kramer, M., & McPhee, C. A. 2006, MNRAS, 373, L50
- Taam & Ricker (2010) Taam, R. E. & Ricker, P. M. 2010, New A, 54, 65
- Taam & Sandquist (2000) Taam, R. E. & Sandquist, E. L. 2000, ARA&A, 38, 113
- Taylor & Weisberg (1982) Taylor, J. H. & Weisberg, J. M. 1982, ApJ, 253, 908
- Thorsett et al. (2005) Thorsett, S. E., Dewey, R. J., & Stairs, I. H. 2005, ApJ, 619, 1036
- Timmes et al. (1996) Timmes, F. X., Woosley, S. E., & Weaver, T. A. 1996, ApJ, 457, 834
- Tutukov & Yungel’Son (1993) Tutukov, A. V. & Yungel’Son, L. R. 1993, Astronomy Reports, 37, 411
- van den Heuvel (2007) van den Heuvel, E. P. J. 2007, in American Institute of Physics Conference Series, Vol. 924, The Multicolored Landscape of Compact Objects and Their Explosive Origins, ed. T. di Salvo, G. L. Israel, L. Piersant, L. Burderi, G. Matt, A. Tornambe, & M. T. Menna, 598–606
- van den Heuvel & Heise (1972) van den Heuvel, E. P. J. & Heise, J. 1972, Nature Physical Science, 239, 67
- Webbink (1984) Webbink, R. F. 1984, ApJ, 277, 355
- Webbink (2008) Webbink, R. F. 2008, in Astrophysics and Space Science Library, Vol. 352, Astrophysics and Space Science Library, ed. E. F. Milone, D. A. Leahy, & D. W. Hobill, 233–+
- Wex et al. (2000) Wex, N., Kalogera, V., & Kramer, M. 2000, ApJ, 528, 401
- Willems & Kalogera (2004) Willems, B. & Kalogera, V. 2004, ApJ, 603, L101
- Willems et al. (2004) Willems, B., Kalogera, V., & Henninger, M. 2004, ApJ, 616, 414
- Willems et al. (2006) Willems, B., Kaplan, J., Fragos, T., Kalogera, V., & Belczynski, K. 2006, Phys. Rev. D, 74, 043003
- Wolszczan (1991) Wolszczan, A. 1991, Nature, 350, 688
- Wong et al. (2010) Wong, T., Willems, B., & Kalogera, V. 2010, ApJ, 721, 1689
- Woods & Ivanova (2011) Woods, T. E. & Ivanova, N. 2011, ApJ, 739, L48
- Woods et al. (2011) Woods, T. E., Ivanova, N., van der Sluys, M., & Chaichenets, S. 2011, in Astronomical Society of the Pacific Conference Series, Vol. 447, Evolution of Compact Binaries, ed. L. Schmidtobreick, M. R. Schreiber, & C. Tappert, 127
|Model||ECS RangeaaThe He core mass range within which a partially degenerate CO core forms. See §2.1.1.||KickbbThe kick velocity applied to a NS born in an Fe core-collapse SN is drawn randomly from a Maxwellian distribution with this dispersion velocity.||HCEccIf on, NS accrete hypercritically in the CE.||CEddIf on, HG stars are allowed to survive a CE.||Notesee: the specific angular momentum of matter; : the maximum mass of a NS; : He stars below this mass develop convective envelopes; F: The fraction of mass lost during stable mass transfer; : CE evolution is determined based on angular momentum conservation, not energy; Tides: Tidal dissipation is decreased by a factor of 5; Twin Binaries: The initial mass function for the secondary is chosen so that the mass ratio is closer to 1; Wind 1: He star winds are decreased by a factor of 4; Wind 2: H and He star winds are decreased by a factor of 4; V: Kick velocities for NS born in ECS are drawn from Maxwellian distributions with this dispersion velocity in km s instead of 1/10th that of Fe core-collapse SN.||log()||log()||RankffThe rank of each model when ordered by . Models with rank 1-20 are in bold.|
|42||1.83-2.5||0.3||300||off||off||F = 1||-23.5||-24.7||59|
|50||1.83-2.5||1.0||300||off||off||F = 1||-23.3||-50.0||118|
|54||1.83-2.5||1.0||300||off||off||F = 0||-26.2||-26.9||81|
|88||2.0-2.5||0.5||300||off||off||F = 0.1||-23.9||-24.5||55|
|89||2.0-2.5||0.5||300||off||off||F = 0||-23.5||-24.3||50|
|90||2.0-2.5||0.5||300||off||off||F = 0.3||-23.8||-24.2||45|
|91||2.0-2.5||0.5||300||off||off||F = 1||-23.8||-50.0||111|
|101||2.0-2.5||0.5||50||off||off||V = 50||-31.7||-33.2||93|
|105||2.0-2.5||1.0||300||off||off||V = 60||-50.0||-50.0||108|
|106||2.0-2.5||1.0||300||off||off||V = 18||-50.0||-50.0||106|
|107||2.0-2.5||1.0||300||off||off||V = 6||-50.0||-50.0||107|
|108||2.0-2.5||1.0||300||off||off||V = 0||-50.0||-50.0||105|
|109||2.0-2.5||1.0||300||on||off||F = 0.1||-24.8||-26.8||78|
|110||2.0-2.5||1.0||300||on||off||F = 1||-23.1||-50.0||104|
|133||2.2-3.0||1.0||300||off||off||F = 0.1||-28.0||-29.2||90|
|134||2.2-3.0||1.0||300||off||off||F = 1||-26.8||-50.0||96|