Although -body simulations have the potential to provide the most detailed population syntheses of relativistic binaries in globular clusters, there are very few actual populations described in the literature. Most of the current work that treats binaries in a consistent and detailed way is limited to open clusters [190, 113, 112, 139, 149] and is focused on a particular outcome of the binary population, such as blue stragglers , brown dwarfs , initial binary distributions , or white dwarf CMD sequences . Portegies Zwart et al. focus on photometric observations of open clusters  and on spectroscopy . In their comparison of -body and Fokker–Planck simulations of the evolution of globular clusters, Takahashi and Portegies Zwart  followed the evolution of = 1 K, 16 K, and 32 K systems with initial mass functions given by Equation (9) and initial density profiles set up from King models. Although they allowed for realistic stellar binary evolution in their comparisons, their focus was on the structural evolution of globular clusters. Consequently there is no binary population provided. Other -body simulations suffer from this same problem . On the other hand, recent work by Shara and Hurley has focused specifically on white dwarf binary populations in globular clusters and has produced a detailed table of close white dwarf binaries that were generated in their simulation .
It is possible to generate a population distribution for black hole binaries in globular clusters using the -body simulations of Portegies Zwart and McMillan  that were intended to describe the population of black hole binaries that were ejected from globular clusters. Their scenario for black hole binary ejection describes a population of massive stars that evolves into black holes. The black holes then rapidly segregate to the core and begin to form binaries. As the black holes are significantly more massive than the other stars, they effectively form a separate sub-system, which interacts solely with itself. The black holes form binaries and then harden through binary-single black hole interactions that occasionally eject either the binary, the single black hole, or both.
They simulated this scenario using = 2048 and = 4096 systems with 1% massive stars. The results of their simulations roughly confirm a theoretical argument based on the recoil velocity that a binary receives during an interaction. Noting that each encounter increases the binding energy by about 20% and that roughly 1/3 of this energy goes into binary recoil, the minimum binding energy of an ejected black hole binary is
At the end of this phase of black hole binary ejection, there is a 50% chance that a binary remains in the cluster with no other black hole to eject it. Thus, there should be a stellar mass black hole binary remaining in about half of the galactic globular clusters. The maximum binding energy of the remaining black hole binary is and is also given by Equation (30). We can then approximate the distribution in energies of the remaining black hole binaries as being flat in . The eccentricities of this population will follow a thermal distribution with .
Dynamical Monte Carlo simulations can be used to study the evolution of binary populations within evolving globular cluster models. Rasio et al.  have used a Monte Carlo approach (described in Joshi et al. [125, 126]) to study the formation and evolution of NS–WD binaries, which may be progenitors of the large population of millisecond pulsars being discovered in globular clusters (see Section 3.3). In addition to producing the appropriate population of binary millisecond pulsars to match observations, the simulations also indicate the existence of a population of NS–WD binaries (see Figure 9).
The tail end of the systems in group B of Figure 9 represents the NS–WD binaries that are in very short period orbits and are undergoing a slow inspiral due to gravitational radiation. These few binaries can be used to infer an order of magnitude estimate on the population of such objects in the galactic globular cluster system. If we consider that there are two binaries with orbital period less than 2000 s out of in 47 Tuc, and assume that this rate is consistent throughout the globular cluster system as a whole, we find a total of 60 such binaries. Although this estimate is quite crude, it compares favorably with estimates arrived at through the encounter rate population syntheses, which are discussed in Section 5.3.3.
More recent applications of the Monte Carlo simulations that have focused on the properties of binaries include Fregeau et al.  who look at the production of blue stragglers and other collision products as a result of binary interactions in globular clusters and Ivanova et al.  who have studied the evolution of binary fractions in globular clusters. The latter work demonstrates the gradual burning of binaries in the core that delays the collapse of the core. In addition, they have also shown the build-up of short period white dwarf binaries in the core through dynamical interactions (see Figure 10).
There is also great promise for the hybrid gas/Monte Carlo method being developed by Spurzem and Giersz . Their recent simulation of the evolution of a cluster of 300,000 equal point-mass stars and 30,000 binaries yields a wealth of detail about the position and energy distribution of binaries in the cluster . Further improvements on their code have resulted in direct integration of the binary-binary and binary-single interactions . As a result, they have been able to produce empirical cross-sections for eccentricity variations during interactions.
One method for exploring the production of relativistic binary populations in globular clusters involves determining the encounter rate expected between different classes of objects in a globular cluster. Sigurdsson and Phinney  use Monte Carlo simulations of binary encounters to infer populations using a static background cluster described by an isotropic King–Michie model. Their results are focused toward predicting the observable end products of binary evolution such as millisecond pulsars, cataclysmic variables, and blue stragglers. Therefore, there are no clear descriptions of relativistic binary populations provided. The work of Ivanova et al.  also uses this technique to determine the evolution of binary fractions, but they also do not provide sufficient detail of the population to distinguish the relativistic binaries from other binaries in the simulation. There is promise to produce a more detailed description of ultracompact X-ray binaries consisting of a white dwarf and a neutron star using encounter rates .
Davies and collaborators use the technique of calculating encounter rates (based on calculations of cross-sections for various binary interactions and number densities of stars using King–Michie static models) to determine the production of end products of binary evolution [45, 43]. Although they also do not provide a clear description of a population of relativistic binaries, their results allow the estimation of such a population.
Using the encounter rates of Davies and collaborators [45, 43], one can follow the evolution of binaries injected into the core of a cluster. A fraction of these binaries will evolve into compact binaries which will then be brought into contact through the emission of gravitational radiation. By following the evolution of these binaries from their emergence from common envelope to contact, we can construct a population and period distribution for present day globular clusters . For a globular cluster with dimensionless central potential , Davies  followed the evolution of 1000 binaries over two runs. The binaries were chosen from a Salpeter IMF with exponent , and the common envelope evolution used an efficiency parameter . One run was terminated after 15 Gyr and the population of relativistic binaries which had been brought into contact through gravitational radiation emission was noted. The second run was allowed to continue until all binaries were either in merged or contact systems. There are four classes of relativistic binaries that are brought into contact by gravitational radiation: low mass WD–WD binaries () with total mass below the Chandrasekhar mass; high mass WD–WD binaries () with total mass above the Chandrasekhar mass; NS–WD binaries (NW); and NS–NS binaries . The number of systems brought into contact at the end of each run is given in Table 3.
In the second run, the relativistic binaries had all been brought into contact. In similar runs, this occurs after another 15 Gyr. An estimate of the present-day period distribution can be made by assuming a constant merger rate over the second 15 Gyr. Consider the total number of binaries that will merge to be described by . Thus, the merger rate is . Assuming that the mergers are driven solely by gravitational radiation, we can relate to the present-day period distribution. We define to be the number of binaries with period less than P, and thus
The merger rate is given by the number of mergers of each binary type per 1000 primordial binaries per 15 Gyr. If the orbits have been circularized (which is quite likely if the binaries have been formed through a common envelope), the evolution of the period due to gravitational radiation losses is given by 
Following this reasoning and using the numbers in Table 3, we can determine the present day population of relativistic binaries per 1000 primordial binaries. To find the population for a typical cluster, we need to determine the primordial binary fraction for globular clusters. Estimates of the binary fraction in globular clusters range from 13% up to about 40% based on observations of either eclipsing binaries [4, 242, 243] or luminosity functions [204, 205]. Assuming a binary fraction of 30%, we can determine the number of relativistic binaries with short orbital period for a typical cluster with and the galactic globular cluster system with  by simply integrating the period distribution from contact up to ,, .
The expected populations for an individual cluster and the galactic cluster system are shown in Table 4 using neutron star masses of , white dwarf masses of and , and .
Although we have assumed the orbits of these binaries will be circularized, there is the possible exception of binaries, which may have a thermal distribution of eccentricities if they have been formed through exchange interactions rather than through a common envelope. In this case, Equations (33) and (34) are no longer valid. An integration over both period and eccentricity, using the formulae of Pierro and Pinto , would be required.
The small number of observed relativistic binaries can be used to infer the population of dark progenitor systems . For example, the low-mass X-ray binary systems are bright enough that we see essentially all of those that are in the galactic globular cluster system. If we assume that the ultracompact ones originate from detached WD–NS systems, then we can estimate the number of progenitor systems by looking at the time spent by the system in both phases. Let be the number of ultracompact LMXBs and be their typical lifetime. Also, let be the number of detached WD–NS systems that will evolve to become LMXBs, and be the time spent during the inspiral due to the emission of gravitational radiation until the companion white dwarf fills its Roche lobe. If the process is stationary, we must have
There are four known ultracompact LMXBs  with orbital periods small enough to require a degenerate white dwarf companion to the neutron star. There are six other LMXBs with unknown orbital periods. Thus, . The lifetime is rather uncertain, depending upon the nature of the mass transfer and the timing when the mass transfer would cease. A standard treatment of mass transfer driven by gravitational radiation alone gives an upper bound of , but other effects such as tidal heating or irradiation may shorten this to [7, 198]. The value of depends critically upon the evolution of the neutron star–main-sequence binary, and is very uncertain. Both and depend upon the masses of the white dwarf secondary and the neutron star primary. For a rough estimate, we take the mass of the secondary to be a typical He white dwarf of mass and the mass of the primary to be . Rather than estimate the typical period of emergence from the common envelope, we arbitrarily choose . We can be certain that all progenitors have emerged from the common envelope by the time the orbital period is this low. The value of can be determined by using Equation (36) and the radius of the white dwarf as determined by Lynden-Bell and O’Dwyer . Adopting the optimistic values of and , and evaluating Equation (38) gives . Thus, we find , which is within an order of magnitude of the numbers found through dynamical simulations (see Section 5.3.2) and encounter rate estimations (see Section 5.3.3).
Current production of ultracompact WD–NS binaries is more likely to arise through collisions of neutron stars with lower mass red giant stars near the current turn-off mass. The result of such a collision is a common envelope that will quickly eject the envelope of the red giant and leave behind the core in an eccentric orbit. The result of the eccentric orbit is to hasten the inspiral of the degenerate core into the neutron star due to gravitational radiation . Consequently, can be significantly shorter . Adopting a value of gives .
Continuing in the spirit of small number statistics, we note that there is one known radio pulsar in a globular cluster NS–NS binary (B2127+11C) and about 50 known radio pulsars in the globular cluster system as a whole (although this number may continue to grow) . We may estimate that NS–NS binaries make up roughly of the total number of neutron stars in the globular cluster system. A lower limit on the number of neutron stars comes from estimates of the total number of active radio pulsars in clusters, giving . Thus, we can estimate the total number of NS–NS binaries to be 2000. Not all of these will be in compact orbits, but we can again estimate the number of systems in compact orbits by assuming that the systems gradually decay through gravitational radiation and thus, and integrating Equation (38) for two neutron stars, we find . Again this value compares favorably with the values found from encounter rate estimations.
This work is licensed under a Creative Commons Attribution-Noncommercial-No Derivative Works 2.0 Germany License.