The formation of very wide binaries during the star cluster dissolution phase
Abstract
Over the past few decades, numerous wide ( au) binaries in the Galactic field and halo have been discovered. Their existence cannot be explained by the process of star formation or by dynamical interactions in the field, and their origin has long been a mystery. We explain the origin of these wide binaries by formation during the dissolution phase of young star clusters: an initially unbound pair of stars may form a binary when their distance in phasespace is small. Using body simulations, we find that the resulting wide binary fraction in the semimajor axis range for individual clusters is , depending on the initial conditions. The existence of numerous wide binaries in the field is consistent with observational evidence that most clusters start out with a large degree of substructure. The wide binary fraction decreases strongly with increasing cluster mass, and the semimajor axis of the newly formed binaries is determined by the initial cluster size. The resulting eccentricity distribution is thermal, and the mass ratio distribution is consistent with gravitationallyfocused random pairing. As a large fraction of the stars form in primordial binaries, we predict that a large number of the observed “wide binaries” are in fact triple or quadruple systems. By integrating over the initial cluster mass distribution, we predict a binary fraction of a few per cent in the semimajor axis range in the Galactic field, which is smaller than the observed wide binary fraction. However, this discrepancy may be solved when we consider a broad range of cluster morphologies.
keywords:
Binaries: general – star clusters – methods: body simulations1 Introduction
A significant fraction of stars in the Galactic field are in binary and multiple systems (e.g, Duquennoy & Mayor, 1991; Fischer & Marcy, 1992; Mason et al., 1998; Shatsky & Tokovinin, 2002; Goodwin & Kroupa, 2005; Kobulnicky & Fryer, 2007; Kouwenhoven et al., 2005, 2007; Lada, 2006; Zinnecker & Yorke, 2007; Goodwin et al., 2007). It is also thought that the majority of stars are born in star clusters (Lada & Lada, 2003). Therefore the majority of binaries^{1}^{1}1For brevity we will use ‘binaries’ to mean ‘multiples’ of any multiplicity for the remainder of this paper, only drawing a distinction where it is necessary. in the field population presumably originate from clustered star formation.
It is well known that binaries are dynamically processed in star clusters with wider and less bound systems tending to be destroyed by encounters (Heggie, 1975; Hills, 1975). Therefore, the field binary population is dynamically processed with respect to the birth population of binaries (Kroupa, 1995; Parker et al., 2009). The origin of most field binaries can be understood as a mixture of differently processed initial populations (Goodwin, 2009).
However, a significant number of very wide ( au) binaries have been observed in the field (see § 2). As such wide binaries are extremely sensitive to destruction they have been used to constrain the properties of the Milky Way. Wide binaries in the Galactic disc and halo have been used to place limits on the density of MACHOs and other unseen material (e.g,. Bahcall et al., 1985; Quinn et al., 2009), to constrain the formation history of the Galaxy (e.g., Allen et al., 2007), and to test the dark matter hypothesis (e.g., Hernandez & Lee, 2008). But the extreme sensitivity to destruction makes the survival (and even formation) of such binaries in a cluster something of a mystery.
Many very wide binaries have separations comparable to the average interstellar separation in clusters (typically a few au), and the very widest binaries have separations of order the size of a young cluster core (typically a few au). Given this, it is difficult to see how they could even form, let alone survive, in a cluster (see, e,g., Scally et al., 1999; Parker et al., 2009). Even in an isolated star forming region the typical size of a star forming core is only au (WardThompson et al., 2007) which presumably sets the very maximum size of a (primordial) binary system.
It is possible that a wide binary forms via dynamical interactions in the Galactic field (the capture mechanism). A prerequisite for this mechanism to work is that a significant amount of kinetic energy is dissipated. This energy dissipation can occur due to tidal friction and due to threebody interaction. Tidal friction occurs in the rare event of two stars nearly colliding in a close encounter. In the vast majority of the cases this results in a merger or a flyby, and only in a small number of cases does this lead to the formation of a binary system. However, all binaries resulting from capture by tidal friction are very tight, with orbital period of several days.
Another possible mechanism is threebody interactions. In this case the third star acts as the energy sink, and is generally ejected with high velocity. However, the stellar density in the field is low, of order , so that threebody encounters are rare, and capture by dynamical friction rarely occurs. Goodman & Hut (1993) find that the creation rate for binaries per unit volume can be approximated by
(1) 
where is the typical mass of a star in the field, the number density of stars, the velocity dispersion, and the gravitational constant. In the solar neighbourhood , and km s. For the field therefore, . This shows that the formation of binaries in the field is extremely rare. Note that in dense star clusters the stellar density is high and the velocity dispersion modest, such that the number of binaries formed via threebody interactions (Eq. 1) may be substantial. However, wide binaries, which have semimajor axes comparable to the size of these clusters, are not formed, as they simply do not fit in these star clusters. Furthermore, body simulations by Kroupa & Burkert (2001) have shown that the observed broad period distribution of binaries in the field cannot be produced by dynamically modifying a tighter period distribution in a star cluster.
A third possibility for the origin of wide binaries is formation during cluster dissolution, which is the mechanism we propose in this paper. In an evolving star cluster, stars that are initially unbound, may become bound to each other as the cluster expands, i.e., if the gravitational influence of the other cluster members decreases^{2}^{2}2Interestingly, Levison et al. (2009) have independently proposed that large populations of comets may be captured by stars during cluster dissolution, by a mechanism which is similar to that discussed in this paper for wide binary formation (see also Eggers et al., 1997).. In order to form a binary pair in this way, (i) the two stars need to be sufficiently close together, (ii) the two stars need to have a sufficiently small velocity difference, and (iii) the newly formed binary should not be destroyed by gravitational interaction with the remaining cluster stars or field stars.
Throughout this paper we refer to binaries with pc as the wide binary population. The paper is organised as follows. In § 2 we briefly discuss surveys of wide binary systems and the corresponding observational techniques. In § 3 we explain our technique and assumptions. In § 4 we provide analytical and Monte Carlo estimates for the resulting wide binary population, and in § 5 we present estimates based on body simulations of evolving star clusters. Finally, in § 6 we present and discuss our conclusions.
2 Observations of wide binaries
Observations have indicated that binaries as wide as 1 pc exist in the halo, while in the Galactic disc the widest binaries have separations of order 0.1 pc (e.g., Close et al., 1990; Chanamé & Gould, 2004). Wider binaries are rare, although some authors claim evidence for binary and higherorder multiple systems wider than 0.1 pc (e.g., Scholz et al., 2008; Caballero, 2009; Mamajek et al., 2009).
Statistical properties of a wide binary population are often recovered using the angular twopoint correlation function (e.g., Bahcall & Soneira, 1981; Garnavich, 1988; Gould et al., 1995; Longhitano & Binggeli, 2010). The most prominent disadvantages of this method are the inability to identity individual wide binaries and the need for a good model for the stellar population that is studied.
Individual wide binary candidates are often identified by their common proper motion on the sky (e.g., Wasserman & Weinberg, 1991; Chanamé & Gould, 2004; Lépine & Bongiorno, 2007; Makarov et al., 2008, and numerous others). Many wide binaries were found by Hipparcos (ESA, 1997) as well; see also Söderhjelm (2007). Their nature can then be further constrained by measuring the parallax and radial velocity (e.g., Latham et al., 1984; Hartigan et al., 1994; Quinn et al., 2009).
The most wellknown survey for binarity is arguably that of Duquennoy & Mayor (1991), who carried out a large binarity study for solartype stars in the solar neighbourhood. They found a lognormal period distribution with a mean and a standard deviation , where is the orbital period in days, in the range days. Note that the latter value roughly corresponds to an orbital period of 30 Myr. In this lognormal period distribution, of the binaries have a semimajor axis wider than au (see Fig. 2).
Several observational studies suggest a semimajor axis distribution of the form , which corresponds to a flat distribution in , also known as Öpik’s law (Öpik, 1924; van Albada, 1968; Vereshchagin et al., 1988; Allen et al., 2000; Poveda & Allen, 2004). In particular, Poveda et al. (2007) find binaries in the field follow Öpik’s law in the separation range au, and suggest that a population of very young binaries follows Öpik’s law up to au (0.2 pc).
Many other authors have also found significant wide binary populations, notably Close et al. (1990), Lépine & Bongiorno (2007), and Chanamé & Gould (2004). We summarise the wide binary separation distribution from various authors in Fig. 1.
The reliability of the observed properties of wide binary populations remains an issue. It is difficult to confidently establish whether the two components of a candidate wide binary are truly bound, or whether it is merely a chance superposition. Due to the long orbital periods of wide binaries, up to millions of years, it is impossible to accurately derive orbital properties and therefore confirm the bound state of the candidate wide binary. It should be noted that due to the lack of detailed orbital information, the quoted separations of wide binaries are usually the instantaneous angular separation. The true separations may thus be significantly larger than the observed, projected separations. On the other hand, an estimate for the semimajor axis distribution of an ensemble of binary systems can statistically be obtained from the projected separation (e.g., Leinert et al., 1993).
On the other hand, it is extremely difficult to identify binaries with distant and faint stellar or substellar companions, due to confusion with foreground and background stars (e.g., Chandrasekhar, 1944). The wide binary fraction as identified in the observational papers may therefore be a lower limit, rather than an upper limit.
Although the exact form of the semimajor axis distribution for very wide binaries is not yet known (see Fig. 1), two features appear to be clear: (i) the binary fraction in the separation range pc is roughly 15%, and (ii) a sharp dropoff in the separation range pc is present, likely due to dynamical destruction of the most weakly bound binary systems.
2.1 Stability of very wide binaries
Whether a wide binary in the Galactic field is stable or not, depends primarily on its semimajor axis . Using a Monte Carlo approach, Weinberg et al. (1987) show that binaries with pc are able to survive in the Galactic disk for Gyr, roughly the age of the Galaxy, but they do not find a sharp cutoff in the semimajor axis distribution. On the other hand, a sharp dropoff is observed for binaries at pc due to interactions with other stars, molecular clouds, and the Galactic tidal field (e.g,. Bahcall & Soneira, 1981; Retterer & King, 1982; Mallada & Fernandez, 2001; Jiang & Tremaine, 2009). Gould & Eastman (2006) explain the slopechange at au in the separation distribution of Lépine & Bongiorno (2007) as the result of dynamical interactions in the Galactic field. Very lowmass binary systems are more weakly bound than their solarmass analogues; their typical separation is therefore expected to be considerably smaller (e.g., Burgasser et al., 2003, see also Kraus & Hillenbrand 2009b), although several very low mass binaries with separations up to pc have been detected (e.g., Radigan et al., 2009, and references therein).
3 Method and assumptions
Model  Structure  (pc)  

P1v  Plummer  1/2  
P1e  Plummer  3/2  
P2v  Plummer  1/2  
P2e  Plummer  3/2  
F1v  Fractal,  1/2  
F1e  Fractal,  3/2  
F2v  Fractal,  1/2  
F2e  Fractal,  3/2 
Given the apparent difficulty of forming extremely wide binaries in star clusters, and the even greater difficulty in keeping them bound we suggest that these extremely wide binaries form during the dissolution of a cluster into the field.
Many young star clusters do not survive for more than a few Myr (Lada & Lada, 2003; Fall et al., 2005; Mengel et al., 2005; Bastian et al., 2005). Their rapid destruction is probably due to the expulsion of the residual gas leftover after star formation, which dramatically changes the cluster potential (Goodwin & Bastian, 2006; Goodwin, 2009, and references therein). Stars that are unbound in a cluster potential may become bound to each other after dissolution as the local density decreases, thus forming an extremely wide binary.
In this section we describe the cluster models used in the remainder of this paper. We use two approaches to investigate the formation of wide binaries during cluster dissolution: Monte Carlo simulations (§ 4.2) and body simulations (§ 5). In § 3.1 we explain our choices for the models used in our analysis, in § 3.2 we define the quantities we use to describe binarity, and in § 3.3 we describe the algorithms we use to ensure the stability of the newly formed wide binaries.
3.1 Model setup
We simulate star clusters using the STARLAB package (Portegies Zwart et al., 2001). We draw single stars from the Kroupa (2001) mass function, , in the mass range . The lowerlimit corresponds to the hydrogenburning limit. The upper limit is (somewhat arbitrarily) set to .
We perform simulations with varying , ranging from small stellar systems (or subclumps) with to open clustersized systems (). We additionally perform simulations of clusters with different radii ( pc), to identify the relation between the initial cluster size and the properties of the newly formed wide binary population. We study two sets of dynamical models: Plummer models and substructured (fractal) models, which we describe below. The properties for the subsets of models are listed in Table 1.
(i) Plummer models. The Plummer sphere is often used in star cluster simulations. In this model, each star is given a certain position and velocity according to the Plummer model (Plummer, 1911) with a certain virial radius . The models are assigned virial radii of pc, which are typical for young clusters. Plummer models are isotropic, and the stellar velocities follow roughly a Maxwellian distribution (Fig. 3).
(ii) Fractal models. Young star clusters show a significant fraction of substructure (e.g., Larson, 1995; Elmegreen, 2000; Testi et al., 2000; Lada & Lada, 2003; Gutermuth et al., 2005; Allen et al., 2007). We set the fractal dimension to (fractal). For comparison, a value corresponds to a homogeneous sphere with radius . Each star is assigned a velocity, as described in Goodwin & Whitworth (2004), such that nearby stars have similar velocities. As in the Plummer models, each cluster is assigned a radius in the range pc. Note, however, the difference between the definition of for the two sets of models.
The virial ratio of a star cluster is defined as the ratio between its kinetic energy and potential energy . Clusters with are in virial equilibrium, and those with and are contracting and expanding, respectively. We study both clusters in virial equilibrium (), as well as clusters with . The latter value for is expected for young clusters with an effective star forming efficiency of 33% (Goodwin & Bastian, 2006). We perform the simulations until the clusters are completely dissolved, which is typically of the order of Myr, the timescale at which the majority of lowmass star clusters are destroyed (see, e.g., Tutukov, 1978; Boutloukos & Lamers, 2003; Bastian et al., 2005; Fall et al., 2005, and numerous others).
3.2 Binarity and multiplicity
At first, we will consider star clusters that initially consist of single stars only, while later (§ 5.3) we will also include primordial binaries. We do not study the evolution of star clusters with primordial higherorder () multiple systems. However, these higherorder systems do form in our star cluster simulations. In this case, the following three useful quantities describing the multiplicity of a stellar population can be used:
(2)  
(3)  
(4) 
(see, e.g., Reipurth & Zinnecker, 1993; Kouwenhoven et al., 2005). Here, , , and denote the number of single stars, binaries, and triples in the system. The quantity is the multiplicity fraction (commonly known as the “binary fraction”). is the nonsingle star fraction, as is the fraction of stars that are single. Finally, is the companion star fraction, which describes the average number of companions per system, where “system” can refer to a single star or multiple system. The number of systems is given by , while denotes the total number of individual stars. Clusters that only contain single stars and binary systems have .
3.3 Detection and stability of wide binary and multiple systems
After each simulation, potential binary and multiple systems are identified as those pairs with negative energy (see also Parker et al., 2009). A multiple () system can only survive for a considerable amount of time if (i) the system is internally stable, and (ii) if the outer orbit is stable against perturbations and tidal forces in the Galactic field. To ensure internal stability of each level in the hierarchy of the multiple system, we impose the Valtonen stability criterion , where and are the semimajor axes of the inner and outer orbits. Valtonen et al. (2008) find that
(5) 
where is the (total) mass of the inner component, the mass of the outer component, the relative inclination of the orbits, and the eccentricity of the outer orbit. For a typical system consisting of equalmass stars, a circular outer orbit () and a prograde outer orbit , the above expression reduces to . Systems with are internally stable for at least revolutions of the outer component. For wide binaries, with orbital periods of years ( au), this corresponds to an internally stable period of at least 2 Gyr.
Wide orbits may additionally be unstable against the tidal forces in the Galactic field and interactions with other single stars and binaries. We therefore additionally impose a maximum semimajor axis of 0.1 pc on the outermost orbit of a binary or multiple, motivated by the observed wide binary population (see § 2). The stability of wider binaries is difficult to assess. As several binaries wider than 0.1 pc are known, our predictions may slightly underestimate the wide binary fraction. The properties of the wide binary populations described in this paper therefore pertain to binaries in the separation range pc. Note that these binary systems fall well in the category “extremely wide binaries” in the Zinnecker (1984) classification of orbital separations.
4 Analytic and numerical estimates
Before proceeding to the body simulations in § 5, it is useful to first obtain some analytical approximations for the prevalence of wide binaries that form during cluster dissolution as well as their orbital characteristics. To this end, we first obtain rough estimates using an analytical approach (§ 4.1), and subsequently using a Monte Carlo approach (§ 4.2).
4.1 Binary formation in Maxwellian velocity space
Wide binaries may form during the dissociation of a cluster if the relative velocity between two stars is sufficiently small that they become bound once the perturbing cluster potential is removed. Here we assume that a wide binary will form if the relative velocity is smaller than or roughly equal to the orbital velocity of a star in a wide binary system. For two stars with masses and in a circular binary orbit with semimajor axis , the velocities of the individual stars are given by
(6) 
where is the gravitational constant and is the mass ratio of the binary system. When adopting, for simplicity, , the above expression reduces to:
(7) 
where is the velocity of either of the components. In order to be able to form a binary with semimajor axis , we require velocity differences to be smaller than the critical velocity
(8) 
For our choice of the IMF, the total mass of a binary system is of order . Binaries with , and au thus typically require velocity differences of , and km s, respectively.
If the velocity distribution of the stars in a given star cluster follows a MaxwellBoltzmann distribution, then so does also the distribution of relative speeds between the stars. We define the relative velocity, , with components and magnitude . The distribution over relative speeds is then given by:
(9) 
(Binney & Tremaine, 1987, p. 485), where is the onedimensional velocity dispersion. In the Plummer model is given by:
(10) 
(Heggie & Hut, 2003) for a cluster in virial equilibrium (i.e. ). Here, is the total mass of the cluster, and the distance to the cluster centre. We can rewrite Eq. (10) in units more suitable for the clusters considered in this paper. First, we set , where is the number of stars in the cluster and is their average mass. Using the Kroupa (2001) IMF, . Evaluating the velocity dispersion at the (intrinsic) halfmass radius, of the cluster, we find:
(11) 
The kinetic energy for a star cluster with is three times that of a cluster with , and therefore the corresponding velocity dispersion is
(12) 
As an example, Fig. 3 shows the distribution of velocities for a Plummer model with stars, a virial radius pc, and . The distribution of relative speeds between random pairs of stars in the cluster is shown in the righthand panel. The latter distribution is well approximated by Eq. (9) with km s, the velocity dispersion at the halfmass radius is given by Eq. (11).
To find the relative fraction of pairs in a given star cluster which has a relative speed such that they may become bound when the cluster disperses, we integrate Eq. (9) between and , where is the critical velocity difference (Eq. 8), below which we assume that two stars may become bound after cluster dissolution:
(13) 
where we normalised the fraction to unity by dividing by the integral of Eq. (9) between 0 and . To find the number of pairs with relative speed less than we multiply by . Hence:
(14) 
If is smaller than unity one might expect that the binary fraction is proportional to . In situations where is larger (i.e., larger than unity) and hence many stars are close to each other in velocity space, we might expect to have some competition between the stars to stay bound.
As an example we show in Fig. 4 how varies with the number of stars, , in clusters with and . We show results in the velocity range km s. Depending on their masses and mass ratios, these velocities correspond to binary systems with semimajor axes from pc down to au. Velocities of km s (not shown in Fig. 4) roughly correspond to binary systems with au.
The values of in Fig. 4 are rather high, as compared to wide binary fractions derived in § 4.2 and 5, mainly because also contains companion stars outside the separation range pc considered throughout this paper. In addition, we believe that the predicted values will drop further due to the inefficiency of the process. For example, it is not likely (but also not impossible) that two stars with nearly the same velocity will form a wide binary system, if there are other stars in between them. Furthermore, we have only considered the relative velocities, while we have ignored the relative positions between the stars. However, this analysis does provide a strong upper limit on the (wide) binary fractions we might expect after cluster dissolution.
4.2 Upper limits from a MonteCarlo approach
In the previous section we obtained rough estimates for the number of newly formed binaries using a Maxwellian velocity distribution. However, we were unable to recover the distributions of orbital properties, such as the semimajor axis distribution, and we were not able to take into account the mass spectrum of stars in the cluster. We therefore use a somewhat more sophisticated Monte Carlo approach, to obtain estimates for these properties as a function of cluster size, structure, number of stars, and virial ratio, for the models listed in Table 1. Estimates of these properties are obtained using an ensemble of initial condition snapshots of each model.
We identify the potential binaries in each star cluster as follows. For each star with mass we determine its nearest neighbour. More precisely, we determine the “most bound” neighbour, i.e., the neighbour with mass for which the internal binding energy
(15) 
is most negative. Here, is the distance between the two stars, their velocity difference, and the gravitational constant. Subsequently, we select those pairs of stars that are each other’s mutual nearest neighbours, and assume that they will form a binary system with a semimajor axis of approximately . Note that not all of these bound pairs may actually form a binary system, as their velocities are perturbed by neighbouring stars. We also ignore the possible presence of triple systems and higherorder systems that may form. The results, shown in Fig. 5, should therefore be considered as a firstorder approximation. We will discuss this figure in detail in § 5, where we will compare the results with those of body simulations (shown in Fig. 6).
Based on a simple Monte Carlo approach, we find that the wide binary fraction decreases with increasing stellar density, and mildly decreases with increasing virial ratio. However, several simplifications have been made, and therefore these results have to be interpreted with care. In particular, we have ignored the interaction of each star with all other stars; we have ignored twobody interactions as well as the tidal field of the cluster. In the following section we perform a more accurate analysis to obtain the abundance and properties of wide binaries formed during cluster dissolution, by performing body simulations. We will discuss all properties in detail, and compare these to the results obtained using the analytical and MonteCarlo approaches.
5 Results from body simulations
The previous two sections have shown that wide binary formation during cluster dissolution may well be possible. In particular, small dense clusters seem the most likely to form wide binaries.
In this section we use body simulations of evolving star clusters to study how the properties of the newly formed binaries and multiple systems depend on the initial properties of the clusters. In § 5.1 we describe the orbital properties and multiplicity fractions of wide binaries resulting from a typical star cluster. In § 5.2 we show how the results depend on the initial size and the number of stars in clusters consisting of initially single stars, for the models listed in Table 1. In § 5.3 we study how the results are affected by the presence of primordial binaries.
5.1 Properties of the newly formed binary population
We perform body simulations of Plummer and fractal clusters consisting of stars with radii pc (the virial radius for the Plummer models, or the total radius of fractal models), with initial virial ratios and (i.e., models P2v, P2e, F2v and F2e in Table 1). Fifty realisations of each model are performed to improve the statistics. The resulting distributions over mass, mass ratio, and semimajor axis for the resulting binary population are shown in Figs. 7 (Plummer models) and 8 (fractal models).
The lefthand panels in Figs. 7 and 8 show the separation distribution of the resulting binary population. Note that binaries of all separations are included in these figures, irrespective of whether they are actually able to survive in the Galactic field or not. The figures illustrate that the separation distribution of the newly formed binaries is bimodal, and consists of a smallseparation dynamical peak and a largeseparation dissolution peak. The tighter binaries in the dynamical peak are formed during dynamical encounters in the cluster, and most of them remain mutually bound during the further evolution of the cluster. The wide binaries in the dissolution peak, on the other hand, are formed during the dissolution phase of star clusters.
The two sets of models with a Plummer density distribution result in a small dynamical peak, indicating that dynamical interactions during the lifetime of the cluster generally do not result in the formation of close binaries. This is not surprising, as all stars in the Plummer models are initially given random velocities. Due to the immediate expansion of the Plummer model with , the dynamical peak is completely absent in this case.
For the models with a fractal density distribution in Fig. 8, there are numerous binaries in the dynamical peak. The fact that the dynamical peak is stronger for models F2v/F2e than for models P2v/P2e is due to both the initial positions and the initial velocities being correlated in the fractal models. Although the average distance between two random stars is similar for both models, the average distance between nearest neighbours in the fractal models is smaller (as they are clumpy). As nearby stars in the clumpy structure also have similar velocities, frequent dynamical interactions occur, resulting in a strong dynamical peak.
For our choice of initial conditions, binaries in the dynamical peak have separations in the range au, with a median value near au. The median value is set by the typical distance between stars in the most densely populated regions of the cluster during the formation of these binary systems. Interestingly, this also corresponds to the observed peak in the TaurusAuriga binary separation distribution (e.g., Leinert et al., 1993; Kroupa et al., 1999).
Binaries in the dissolution peak have a semimajor axis in the separation range pc. The widest binaries in the dissolution peak will immediately break up in the Galactic field, hence our choice to study the wide binary population in the separation range pc throughout this paper. The median separation of binaries in the dissolution peak occurs at pc. As we will see later (§ 5.2.2), this value is set by the initial size of the cluster.
Model  

Separation range  au  pc  pc 
P2v ()  0.2%  0.3%  0.8% 
P2e ()  0.1%  1.4%  2.8% 
F2v ()  3.3%  1.8%  2.6% 
F2e ()  2.2%  0.6%  1.1% 
For practical purposes, we consider three ranges in semimajor axis: close binaries with au, wide binaries with pc, and extremely wide binaries with pc. The limits are indicated with the vertical dotted lines in the figures. Most close binaries that are found in star clusters are formed via the “normal” star formation process, with the small number seen in these simulations formed by dynamical interactions. The wide and extremely wide binaries are formed during the cluster dissolution phase. Note however, that the vast majority of the extremely wide binaries are unstable in the Galactic field, and are ionised quickly after their formation.
For the models in Figs. 7 and 8, the specific binary fraction (i.e., the fraction of binary systems in a certain semimajor axis range) of the three types of binaries are listed in Table 2. The highest wide binary fractions of a few per cent (in the separation range pc) are obtained for Plummer models with , and fractal models with .
The middle and righthand panels of Figs. 7 and 8 show the correlations between semimajor axis, mass ratio, and primary mass, for the binary and multiple (higherorder) systems in each of the models. The panels indicate the presence of a large number of newly formed multiple systems. These higherorder systems are stable in isolation, but a large fraction will not be able to survive in the Galactic field, where tidal forces will rapidly remove the outer component from the system. Figs. 7 and 8 therefore overestimate the fraction of higherorder multiple systems. Note, in particular, the high prevalence of multiple systems in the dynamical peak. Many outer components of these systems fall in the dissolution peak. These systems are therefore wide higherorder systems.
For the Plummer models, the dynamical peak consists of systems with high masses and high mass ratios. This is a wellknown signature of mass segregation: the highestmass stars sink to the cluster centre, where they form close binaries (e.g. Heggie & Hut, 2003). During the dissolution phase of the clusters, these close, massive binaries act like single stars when forming a “wide binary”, which is in fact a wide triple or higherorder multiple system. The effect of mass segregation is less visible for the fractal models, where dynamical interactions in the subclumps play a greater role. However, Fig. 8 still clearly shows that most massive systems are mostly close ( au) and often higherorder. In addition to the triple and higherorder systems formed during the dissolution process, several higherorder systems may form via dynamical interactions (van den Berk et al., 2007).
The correlations between the primary mass and mass ratio distributions for the binary and multiple systems are similar to those expected from random pairing of individual stars. To first order approximation, the masses of the two stars, and , in each binary are uncorrelated; one would therefore expect something similar to random pairing (e.g., Kouwenhoven et al., 2009), where the average mass ratio decreases with increasing binary system mass. The resulting mass ratio distributions for binary and multiple systems with au (i.e., those in the dissolution peak) are shown in Fig. 9, which illustrates the dependence of the mass ratio distribution on mass.
Based on the analysis of a sample of 798 common proper motion pairs, Trimble (1987) also come to the conclusion that the very wide binary population in the field is consistent with random pairing, and Valtonen (1997) come to the same conclusion from their simulations of threebody encounters. However, the wide binary population does not result from random pairing alone, as the interaction between two stars depends on their mutual gravitational attraction, and the probability of two stars forming a binary is thus proportional to the product . In other words, gravitational focusing (e.g., Gaburov et al., 2008) plays an important role.
Measurements of the eccentricity distribution of wide binaries are currently unavailable, due to the large orbital periods and incompleteness. If we suspect that the vast majority of wide binaries probably have formed dynamically, and as dynamical interactions are common among the widest binaries (with respect to closerin binaries), the best guess is perhaps the thermal eccentricity distribution (Heggie, 1975, see also Kroupa 2008 for a derivation), which results from energy equipartition. The eccentricity distributions resulting for binaries in the dissolution peak ( au) are shown in Fig. 10. As expected, the thermal eccentricity distribution is a good approximation for the newly formed binary population.
If wide binaries form during the dissolution process of a star cluster, then the orbital and spin angular momenta of the components should be randomly aligned. On the other hand, if the two components formed together in some way it might be that the orbital and spin angular momenta of the components will be correlated (as seen for example in the observations of au Ae/Be binaries by Baines et al. 2006). Therefore, observations of the relative alignments of orbital and spin angular momenta could provide constraints on the possible formation mechanisms of very wide binaries.
Finally, the age difference (between primary and companion star) for a population of wide binaries could provide a clue to their origin (see, e.g., Kraus & Hillenbrand, 2009a). For a star cluster with a certain age spread, one might expect the components of the resulting wide binary population to exhibit a similar age difference. This age difference is measurable, but only for young ( Myr) binary systems. On the other hand, this age difference may be smaller than expected from random pairing, if an initial correlation between position and velocity exists.
5.2 Dependence on cluster properties
In this section we describe how the properties of the wide binary population depend on the initial conditions we assign to a star cluster, in particular its size , number of stars , virial ratio , and morphology (Plummer sphere or fractal structure). We adopt the cluster properties listed in Table 1. We compare the results that we derived earlier using Monte Carlo simulations (Fig. 5), with the results of body simulations, shown in Fig. 6.
5.2.1 Dependence on the initial cluster mass
The top panels of Fig. 6 show the median semimajor axis and the binary fraction of wide binaries ( pc) as a function of the number of stars in a cluster. For both the Plummer models and the fractal models, does not vary significantly with and the virial ratio . The reason for this is that all these models have an identical size . Since is the most important size scale imposed on the modelled star clusters, it determines the size scaling (i.e., semimajor axis distribution) of the newly formed binaries.
The dependence of on and is qualitatively the same as the analytical predictions shown in Fig. 4 and the Monte Carlo approximation shown in Fig. 5 . The wide binary fraction decreases with increasing because the stars are further apart in velocity space (cf. Fig. 3), i.e., the velocity dispersion is larger, and hence two neighbouring stars are less likely to form a bound system.
For the body simulations we find that the fractal model with provides the highest wide binary fractions, although the difference between models is fairly small (especially when compared to the difference with increasing ). Models with generally result in a smaller wide binary fraction than those with , due to the larger distance between the stars in velocity space (see Eqs. 11 and 12). The curves for the fractal models in Figs. 5 and Fig. 6 are almost the same, indicating that the Monte Carlo approach provides a good estimate of the wide binary population. For the Plummer models, the Monte Carlo approach predicts a binary fraction that is too high, which is due to the fact that the positions and velocities of stars in the Plummer models as we initialise them are uncorrelated.
5.2.2 Dependence on the initial cluster size
The bottom panels of Fig 6 shows the dependence of and on the initial size of the clusters. Again, these values are only for wide binaries with pc. Note the different definitions of : for the Plummer models represents the virial radius, while for the fractal models represents the radius of the sphere enclosing the whole system. Note again the similarity between the Monte Carlo approximation shown in Fig. 5 and the body models.
As discussed above, the initial cluster size determines the length scale in each model, and therefore the size scaling of the semimajor axis distribution of the newly formed binaries. For example, changing the initial size of the clusters shown in Figs. 7 and 8 would simply result in the semimajor axis distribution in the lefthand panels being shifted to smaller or larger values of .
This direct dependence of on is not seen directly in Fig. 5 because we only show the results for wide binaries in the separation range pc, and because is bimodal. However, the dependent median semimajor axis and binary fraction can be explained by the dynamical peak and dissolution peak shifting through the range pc whilst varying .
The highest is found when either the dynamical peak, or the dissolution peak, is centred in the separation range pc. For our choice of the initial conditions, this peak occurs at pc for the Plummer models, when the dissolution peak is centred in the range pc. The peak in occurs at pc for the fractal models with and at pc for fractal models with , when the dynamical peak is centred in the range pc.
Given our set of initial conditions, compact clusters result in a wide binary fraction of , irrespective of virial ratio and morphology. For more extended clusters, those with a Plummer structure and those with a higher virial ratio result in a smaller binary fraction. The difference between the Plummer and fractal models can be explained by (i) the difference in the definition of for the two sets of models, and (ii) by the different intrinsic separation distribution (see the lefthand panels in Figs. 7 and 8).
The cluster size determines the length scale of the system, and therefore determines the typical semimajor axis of the newly formed wide binaries. Other, less important length scales in the system are the mean distance between two stars, which depends on the parameters , , and the stellar density distribution (see § 5.2.1), as well as the typical semimajor axis of primordial binary systems (see § 5.3).
5.3 Effects of primordial binarity
In the analysis above we have considered star clusters that initially consist of single stars only. The results for star clusters with a nonzero primordial binary fraction are very similar to the results described above, with the difference that the components of the wide “binary” are now in many cases primordial binaries. In other words, the majority of the wide “binaries” that formed in the simulations described in the previous sections, actually describe the outer orbits of wide triple and quadruple systems.
We predict the properties of these wide triple and quadruple systems by performing body simulations of the Plummer models listed in Table 1, but now we include a nonzero primordial binary fraction. We perform the simulations with primordial binary fractions ranging from 0% to 100%. We adopt the Kroupa (1995) birth period distribution. This distribution is derived from a detailed analysis of observed stellar populations, and has the form
(16) 
for , where , , and is the period in days. We adopt a thermal eccentricity distribution (). We adopt a flat mass ratio distribution with (i.e., we apply pairing function PCPI; see Kouwenhoven et al., 2009). Subsequently, we generate an initial population from this birth population, by applying eigenevolution as described in Kroupa (1995). All binaries are assigned random orientations and orbital phases at the beginning of the simulations.
Due to the inclusion of binary components, the total mass of each cluster increases slightly (up to a maximum of % for a primordial binary fraction of 100%), although the number of “systems”, remains constant. Strictly speaking, it is thus not appropriate to directly compare clusters with and without binaries, as we have changed more than one parameter: binarity and cluster mass (see, e.g., Kouwenhoven & de Grijs, 2008). However, as the increase in cluster mass due to adding the companions is rather small, we will ignore this issue.
The results for clusters with a varying primordial binary fraction is shown in the top panels of Fig. 11. For small binary fractions, the results are very similar to those of clusters without primordial binaries. The properties of the resulting wide binary population depend mildly on . An increasing results in, on average, wider binaries, hence in a larger fraction of binaries with pc, and therefore in a slightly smaller wide binary fraction. Note that decreases slightly with increasing . A larger primordial binary fraction results in a smaller wide binary fraction, possibly because of the destruction of newly formed wide binaries by primordial binaries (which have a significantly larger collisional crosssection than single stars).
Whether or not primordial binarity affects the formation of wide binaries depends not only on the primordial binary fraction, but also on the properties of these binaries: the semimajor axis (or period) distribution, the eccentricity distribution and the mass ratio distribution. The most important of these is the semimajor axis distribution , as it determines the internal binding energy of a binary and the crosssection for gravitational interactions between binaries and other binaries or single stars. In order to extract the dependence on , we simulate clusters in which all binaries have a single value for . We vary in each cluster, and determine the number of newlyformed binaries. In all simulations we adopt a primordial binary fraction of 50%, a flat mass ratio distribution and a thermal eccentricity distribution.
The results of these simulations are shown in the bottom panels of Fig. 11. For models with au, most primordial binaries survive, while additional wide binary, triple, and quadruple systems are formed. In fact, the resulting wide binary fraction is practically independent of the primordial binary fraction . For models with au, all primordial binary systems are classified as wide binaries. For these models we therefore have and a median semimajor axis equal to , which results in the glitches at au in Fig. 11.
The wide orbits are part of systems with , , and components. They are formed by randomly pairing single stars and primordial binary systems together. The number of multiple systems of each degree can thus be estimated by simply calculating the probability of randomly drawing a singlesingle, singlebinary, and binarybinary pair. When assuming a primordial binary fraction , the multiplicity distribution of the resulting wide population can be estimated as follows:
Wide binary fraction  (17)  
Wide triple fraction  (18)  
Wide quadruple fraction  (19) 
where we have made the assumption that none of the primordial binary systems has broken up.
All models shown in Fig. 11 result in wide binary fractions that are more or less independent of and . The value of is therefore primarily determined by the initial values of the number of system in the cluster, and its initial size .
If we assume that the wide orbits in the bottom panels of Fig. 11 (where ) are formed of randomly paired components (i.e., single stars or primordial binaries), we can calculate the fraction of higherorder multiple systems among the wide binaries. Among these, we predict that , , and , are binary, triple, and quadruple systems, respectively. In this example, we thus expect 75% of the “wide binaries” to be higherorder multiple systems. Due to the random process, the outer orbits of these systems are expected to be uncorrelated with the inner orbits or stellar spin axes.
The ratios between wide binary, triple, and quadruple systems are therefore indicative of . A survey for higherorder multiplicity among “wide binary systems” can thus be used to constrain the primordial binary fraction. Given the fact that the majority of stars do form in binary systems, we predict a very high fraction of higherorder multiple systems among wide “binary” systems; see Fig. 12. Our proposed mechanism could explain the existence of the observed wide multiple systems (Mamajek et al., 2009), and our predictions are strongly supported by the surveys of Makarov et al. (2008) and Faherty et al. (2010), who find that a high fraction of the common proper motion pairs in their survey contain inner binary or triple systems, which is significantly higher than in populations of other types of binaries.
6 Summary and discussion
Observations have shown that of binaries are wide ( au). These wide binaries are difficult to explain as being the result of star formation as it is difficult to see how wide binaries can form (especially those au), and they would be rapidly destroyed in clustered star forming environments. Whilst – of stars do appear to form in an ‘isolated’ environment in which such binaries could possibly survive, in order to explain the fraction of wide binaries, almost all stars forming in isolated environments would have to form wide binaries. Further, such wide binaries cannot be formed later in any significant numbers by dynamical interactions in the Galactic field.
In this paper we study the possibility of wide binary formation during the dissolution phase of star clusters, in particular, during the rapid expansion of clusters after gas expulsion. We study this possibility using (1) an analytical approach in an idealised situation, (2) a Monte Carlo approach, and (3) detailed body simulations. Our main conclusions are as follows:

The wide binary fraction among the dissolved stellar population ranges between 1% and 30%, depending on the cluster properties.

More massive star clusters result in a smaller wide than lowmass clusters. Clusters with a spherical, smooth stellar density distribution form fewer wide binaries than substructured clusters of the same size and mass. This is due to the fact that the average distance between nearest neighbours is smaller for substructured clusters. Expanding (postgas expulsion) star clusters produce a larger than those starting out of equilibrium.

The typical semimajor axis of the newly formed binaries is similar to the initial size of the star cluster from which they were born. The resulting semimajor axis distribution is generally bimodal, consisting of a dynamical peak with binary systems formed by dynamical interactions, and a dissolution peak with binary systems formed during the cluster dissolution phase.

The formation of wide binaries during the star cluster dissolution phase is a random process, resulting in the following orbital properties. The eccentricity distribution of the wide binaries is approximately thermal: for . The mass ratio distribution of the wide binaries is the result of gravitationallyfocused random pairing. In a wide binary, the orbital and spin angular momenta are uncorrelated.

Star clusters with a nonzero primordial binary population form wide triple and quadruple systems, i.e., the components of a newlyformed wide “binary” can themselves be close primordial binaries, rather than single stars. The ratio of triple to quadruple systems among very wide orbits is therefore indicative of the primordial binary fraction . Given that is large, we predict a high frequency of triple and quadruple systems among the known wide “binary” systems, which is supported by existing surveys for higherorder multiplicity among wide binary systems.
Throughout this paper we have made predictions of the properties of the wide binary population resulting from the dissolution of individual clusters. In order to compare our results with observations, we should therefore take into account the fact that the field star population is made up of the stars resulting from an ensemble of clusters of different sizes and masses. The initial cluster mass distribution may be approximated by with (see, e.g., Zhang & Fall, 1999; Ashman & Zepf, 2001; Bik et al., 2003; Hunter et al., 2003). Given the number of stars , where is the average mass of a star, this distribution is equivalent to . Oey et al. (2004) suggest that the above expression can be extrapolated down to . The upper limit for the initial cluster mass distribution is (e.g., de Grijs & Parmentier, 2007, and references therein). The resulting binary fraction for the ensemble of stars (i.e., the field star population) is then given by:
(20) 
where is the cluster mass dependent wide binary fraction. The numerator in the above expression is proportional to the number of binaries, and the denominator is proportional to the total number of stars in the ensemble of clusters. In addition, the size and dissolution time of a star cluster, and therefore the wide binary fraction, may also depend on its Galactic location (e.g., Baumgardt & Makino, 2003). An inspection of Fig. 6 shows that an extrapolation to results in a wide binary fraction of several per cent; smaller than the observed 15%, irrespective of the choices for , , and the morphology of the cluster. Although we predict rather small values, our backoftheenvelope calculation does result in the right order of magnitude for the wide binary fraction in the Galactic field. It is clear, however, that a deeper investigation is required to accurately recover the properties of the wide binary population in the field. In particular, a wider range of star cluster morphologies has to be considered, by varying the fractal dimension and positionvelocity correlations of individual star clusters.
Our proposed formation mechanism for very wide binaries predicts at least several common proper motion pairs in and around dissolving star clusters and moving groups. For example, the mechanism may well explain the presence of the three common proper motion pairs in the moving groups studied by Clarke et al. (2009). The future prospects in wide binary research are bright: an enormous number of wide binaries are expected to be found with the GAIA mission^{3}^{3}3http://www.esa.int/science/gaia (Perryman et al., 2001; Turon et al., 2005) and LAMOST ^{4}^{4}4http://www.lamost.org (Chu, 1998; Stone, 2008). These datasets should help determine the true fraction of wide binaries and their orbital parameters.
Acknowledgements
We would like to thank Anthony Whitworth (our referee) and Simon Portegies Zwart for a useful discussion on this topic. MBNK was supported by the Peter and Patricia Gruber Foundation through the PPGF fellowship, the Peking University One Hundred Talent Fund (985), and by PPARC/STFC (grant PP/D002036/1). We would like to acknowledge CICS for the provision of research computing facilities through the White Rose Grid. We acknowledge research support from and hospitality at the International Space Science Institute in Bern (Switzerland), as part of an International Team programme. The authors acknowledge the SheffieldBonn Royal Society International Joint Project grant, which provided financial support and the collaborative opportunities for this work. The calculations performed by DM were carried out on computer hardware purchased from the Royal Physiographic Society in Lund. RJP acknowledges financial support from STFC.
References
 Allen et al. (2007) Allen C., Poveda A., HernándezAlcántara A., 2007, in Hartkopf W. I., Guinan E. F., Harmanec P., eds, IAU Symposium Vol. 240 of IAU Symposium, Halo Wide Binaries and Moving Clusters as Probes of the Dynamical and Merger History of our Galaxy. pp 405–413
 Allen et al. (2000) Allen C., Poveda A., Herrera M. A., 2000, A&A, 356, 529
 Allen et al. (2007) Allen L., Megeath S. T., Gutermuth R., Myers P. C., Wolk S., Adams F. C., Muzerolle J., Young E., Pipher J. L., 2007, in Reipurth B., Jewitt D., Keil K., eds, Protostars and Planets V The Structure and Evolution of Young Stellar Clusters. pp 361–376
 Ashman & Zepf (2001) Ashman K. M., Zepf S. E., 2001, AJ, 122, 1888
 Bahcall et al. (1985) Bahcall J. N., Hut P., Tremaine S., 1985, ApJ, 290, 15
 Bahcall & Soneira (1981) Bahcall J. N., Soneira R. M., 1981, ApJ, 246, 122
 Baines et al. (2006) Baines D., Oudmaijer R. D., Porter J. M., Pozzo M., 2006, MNRAS, 367, 737
 Bastian et al. (2005) Bastian N., Gieles M., Lamers H. J. G. L. M., Scheepmaker R. A., de Grijs R., 2005, A&A, 431, 905
 Baumgardt & Makino (2003) Baumgardt H., Makino J., 2003, MNRAS, 340, 227
 Bik et al. (2003) Bik A., Lamers H. J. G. L. M., Bastian N., Panagia N., Romaniello M., 2003, A&A, 397, 473
 Binney & Tremaine (1987) Binney J., Tremaine S., 1987, Galactic dynamics. Princeton, NJ, Princeton University Press, 1987
 Boutloukos & Lamers (2003) Boutloukos S. G., Lamers H. J. G. L. M., 2003, MNRAS, 338, 717
 Burgasser et al. (2003) Burgasser A. J., Kirkpatrick J. D., Reid I. N., Brown M. E., Miskey C. L., Gizis J. E., 2003, ApJ, 586, 512
 Caballero (2009) Caballero J. A., 2009, A&A, 507, 251
 Chanamé & Gould (2004) Chanamé J., Gould A., 2004, ApJ, 601, 289
 Chandrasekhar (1944) Chandrasekhar S., 1944, ApJ, 99, 54
 Chu (1998) Chu Y., 1998, Highlights of Astronomy, 11, 493
 Clarke et al. (2009) Clarke J. R. A., Pinfield D. J., GálvezOrtiz M. C., Jenkins J. S., Burningham B., Deacon N. R., Jones H. R. A., Pokorny R. S., Barnes J. R., DayJones A. C., 2009, MNRAS, p. 1923
 Close et al. (1990) Close L. M., Richer H. B., Crabtree D. R., 1990, AJ, 100, 1968
 de Grijs & Parmentier (2007) de Grijs R., Parmentier G., 2007, Chinese Journal of Astronomy and Astrophysics, 7, 155
 Duquennoy & Mayor (1991) Duquennoy A., Mayor M., 1991, A&A, 248, 485
 Eggers et al. (1997) Eggers S., Keller H. U., Kroupa P., Markiewicz W. J., 1997, Planet. Space Sci., 45, 1099
 Elmegreen (2000) Elmegreen B. G., 2000, ApJ, 530, 277
 ESA (1997) ESA 1997, VizieR Online Data Catalog, 1239
 Faherty et al. (2010) Faherty J. K., Burgasser A. J., West A. A., Bochanski J. J., Cruz K. L., Shara M. M., Walter F. M., 2010, AJ, 139, 176
 Fall et al. (2005) Fall S. M., Chandar R., Whitmore B. C., 2005, ApJ, 631, L133
 Fischer & Marcy (1992) Fischer D. A., Marcy G. W., 1992, ApJ, 396, 178
 Gaburov et al. (2008) Gaburov E., Gualandris A., Portegies Zwart S., 2008, MNRAS, 384, 376
 Garnavich (1988) Garnavich P. M., 1988, ApJ, 335, L47
 Goodman & Hut (1993) Goodman J., Hut P., 1993, ApJ, 403, 271
 Goodwin (2009) Goodwin S. P., 2009, ArXiv:0911.0795
 Goodwin & Bastian (2006) Goodwin S. P., Bastian N., 2006, MNRAS, 373, 752
 Goodwin & Kroupa (2005) Goodwin S. P., Kroupa P., 2005, A&A, 439, 565
 Goodwin et al. (2007) Goodwin S. P., Kroupa P., Goodman A., Burkert A., 2007, in Reipurth B., Jewitt D., Keil K., eds, Protostars and Planets V The Fragmentation of Cores and the Initial Binary Population. pp 133–147
 Goodwin & Whitworth (2004) Goodwin S. P., Whitworth A. P., 2004, A&A, 413, 929
 Gould et al. (1995) Gould A., Bahcall J. N., Maoz D., Yanny B., 1995, ApJ, 441, 200
 Gould & Eastman (2006) Gould A., Eastman J., 2006, astroph/0610799
 Gutermuth et al. (2005) Gutermuth R. A., Megeath S. T., Pipher J. L., Williams J. P., Allen L. E., Myers P. C., Raines S. N., 2005, ApJ, 632, 397
 Hartigan et al. (1994) Hartigan P., Strom K. M., Strom S. E., 1994, ApJ, 427, 961
 Heggie & Hut (2003) Heggie D., Hut P., 2003, The Gravitational MillionBody Problem: A Multidisciplinary Approach to Star Cluster Dynamics. Cambridge University Press, 2003
 Heggie (1975) Heggie D. C., 1975, MNRAS, 173, 729
 Hernandez & Lee (2008) Hernandez X., Lee W. H., 2008, MNRAS, 387, 1727
 Hills (1975) Hills J. G., 1975, AJ, 80, 809
 Hunter et al. (2003) Hunter D. A., Elmegreen B. G., Dupuy T. J., Mortonson M., 2003, AJ, 126, 1836
 Jiang & Tremaine (2009) Jiang Y., Tremaine S., 2009, MNRAS, p. 1640
 Kobulnicky & Fryer (2007) Kobulnicky H. A., Fryer C. L., 2007, ApJ, 670, 747
 Kouwenhoven et al. (2009) Kouwenhoven M. B. N., Brown A. G. A., Goodwin S. P., Portegies Zwart S. F., Kaper L., 2009, A&A, 493, 979
 Kouwenhoven et al. (2007) Kouwenhoven M. B. N., Brown A. G. A., Portegies Zwart S. F., Kaper L., 2007, A&A, 474, 77
 Kouwenhoven et al. (2005) Kouwenhoven M. B. N., Brown A. G. A., Zinnecker H., Kaper L., Portegies Zwart S. F., 2005, A&A, 430, 137
 Kouwenhoven & de Grijs (2008) Kouwenhoven M. B. N., de Grijs R., 2008, A&A, 480, 103
 Kraus & Hillenbrand (2009a) Kraus A. L., Hillenbrand L. A., 2009a, ApJ, 704, 531
 Kraus & Hillenbrand (2009b) Kraus A. L., Hillenbrand L. A., 2009b, ApJ, 703, 1511
 Kroupa (1995) Kroupa P., 1995, MNRAS, 277, 1507
 Kroupa (2001) Kroupa P., 2001, MNRAS, 322, 231
 Kroupa (2008) Kroupa P., 2008, in S. J. Aarseth, C. A. Tout, & R. A. Mardling ed., Lecture Notes in Physics, Berlin Springer Verlag Vol. 760 of Lecture Notes in Physics, Berlin Springer Verlag, Initial Conditions for Star Clusters. p. 181
 Kroupa & Burkert (2001) Kroupa P., Burkert A., 2001, ApJ, 555, 945
 Kroupa et al. (1999) Kroupa P., Petr M. G., McCaughrean M. J., 1999, New Astronomy, 4, 495
 Lada (2006) Lada C. J., 2006, ApJ, 640, L63
 Lada & Lada (2003) Lada C. J., Lada E. A., 2003, ARA&A, 41, 57
 Larson (1995) Larson R. B., 1995, MNRAS, 272, 213
 Latham et al. (1984) Latham D. W., Schechter P., Tonry J., Bahcall J. N., Soneira R. M., 1984, ApJ, 281, L41
 Leinert et al. (1993) Leinert C., Zinnecker H., Weitzel N., Christou J., Ridgway S. T., Jameson R., Haas M., Lenzen R., 1993, A&A, 278, 129
 Lépine & Bongiorno (2007) Lépine S., Bongiorno B., 2007, AJ, 133, 889
 Longhitano & Binggeli (2010) Longhitano M., Binggeli B., 2010, ArXiv:1001.3819
 Makarov et al. (2008) Makarov V. V., Zacharias N., Hennessy G. S., 2008, ApJ, 687, 566
 Mallada & Fernandez (2001) Mallada E., Fernandez J. A., 2001, in Revista Mexicana de Astronomia y Astrofisica Conference Series Vol. 11 of RMAA Conference Series, Dynamical Evolution of Wide Binaries. p. 27
 Mamajek et al. (2009) Mamajek E. E., Kenworthy M. A., Hinz P. M., Meyer M. R., 2009, ArXiv:0911.5028
 Mason et al. (1998) Mason B. D., Gies D. R., Hartkopf W. I., Bagnuolo Jr. W. G., ten Brummelaar T., McAlister H. A., 1998, AJ, 115, 821
 Mengel et al. (2005) Mengel S., Lehnert M. D., Thatte N., Genzel R., 2005, A&A, 443, 41
 Oey et al. (2004) Oey M. S., King N. L., Parker J. W., 2004, AJ, 127, 1632
 Öpik (1924) Öpik E., 1924, Tartu Obs. Publ., 25, No. 6
 Parker et al. (2009) Parker R. J., Goodwin S. P., Kroupa P., Kouwenhoven M. B. N., 2009, MNRAS, 397, 1577
 Perryman et al. (2001) Perryman M. A. C., de Boer K. S., Gilmore G., Høg E., Lattanzi M. G., Lindegren L., Luri X., Mignard F., Pace O., de Zeeuw P. T., 2001, A&A, 369, 339
 Plummer (1911) Plummer H. C., 1911, MNRAS, 71, 460
 Portegies Zwart et al. (2001) Portegies Zwart S. F., McMillan S. L. W., Hut P., Makino J., 2001, MNRAS, 321, 199
 Poveda & Allen (2004) Poveda A., Allen C., 2004, in Allen C., Scarfe C., eds, Rev. Mexicana Astron. Astrofis.Conference Series The distribution of separations of wide binaries of different ages. p. 49
 Poveda et al. (2007) Poveda A., Allen C., HernándezAlcántara A., 2007, in IAU Symposium Vol. 240 of IAU Symposium, The Frequency Distribution of Semimajor Axes of Wide Binaries: Cosmogony and Dynamical Evolution. pp 417–425
 Quinn et al. (2009) Quinn D. P., Wilkinson M. I., Irwin M. J., Marshall J., Koch A., Belokurov V., 2009, MNRAS, 396, L11
 Radigan et al. (2009) Radigan J., Lafrenière D., Jayawardhana R., Doyon R., 2009, ApJ, 698, 405
 Reipurth & Zinnecker (1993) Reipurth B., Zinnecker H., 1993, A&A, 278, 81
 Retterer & King (1982) Retterer J. M., King I. R., 1982, ApJ, 254, 214
 Scally et al. (1999) Scally A., Clarke C., McCaughrean M. J., 1999, MNRAS, 306, 253
 Scholz et al. (2008) Scholz R.D., Kharchenko N. V., Lodieu N., McCaughrean M. J., 2008, A&A, 487, 595
 Shatsky & Tokovinin (2002) Shatsky N., Tokovinin A., 2002, A&A, 382, 92
 Söderhjelm (2007) Söderhjelm S., 2007, A&A, 463, 683
 Stone (2008) Stone R., 2008, Science, 320, 34
 Testi et al. (2000) Testi L., Sargent A. I., Olmi L., Onello J. S., 2000, ApJ, 540, L53
 Trimble (1987) Trimble V., 1987, Astronomische Nachrichten, 308, 343
 Turon et al. (2005) Turon C., O’Flaherty K. S., Perryman M. A. C., eds, 2005, The ThreeDimensional Universe with Gaia Vol. 576 of ESA Special Publication
 Tutukov (1978) Tutukov A. V., 1978, A&A, 70, 57
 Valtonen et al. (2008) Valtonen M., Mylläri A., Orlov V., Rubinov A., 2008, in IAU Symposium Vol. 246 of IAU Symposium, The Problem of Three Stars: Stability Limit. pp 209–217
 Valtonen (1997) Valtonen M. J., 1997, ApJ, 485, 785
 van Albada (1968) van Albada T. S., 1968, Bull. Astron. Inst. Netherlands, 20, 47
 van den Berk et al. (2007) van den Berk J., Portegies Zwart S. F., McMillan S. L. W., 2007, MNRAS, 379, 111
 Vereshchagin et al. (1988) Vereshchagin S., Tutukov A., Iungelson L., Kraicheva Z., Popova E., 1988, Ap&SS, 142, 245
 WardThompson et al. (2007) WardThompson D., André P., Crutcher R., Johnstone D., Onishi T., Wilson C., 2007, in B. Reipurth, D. Jewitt, & K. Keil ed., Protostars and Planets V An Observational Perspective of LowMass Dense Cores II: Evolution Toward the Initial Mass Function. pp 33–46
 Wasserman & Weinberg (1991) Wasserman I., Weinberg M. D., 1991, ApJ, 382, 149
 Weinberg et al. (1987) Weinberg M. D., Shapiro S. L., Wasserman I., 1987, ApJ, 312, 367
 Zhang & Fall (1999) Zhang Q., Fall S. M., 1999, ApJ, 527, L81
 Zinnecker (1984) Zinnecker H., 1984, Ap&SS, 99, 41
 Zinnecker & Yorke (2007) Zinnecker H., Yorke H. W., 2007, ARA&A, 45, 481