On the properties of free floating planets originating in circumbinary planetary systems

Free-floating planets are a new class of planets recently discovered. These planets don't orbit within stellar systems, instead living a nomadic life within the galaxy. How such objects formed remains elusive. Numerous works have explored mechanisms to form such objects, but have not yet provided predictions on their distributions that could differentiate between formation mechanisms. In this work we form these objects within circumbinary systems, where these planets are readily formed and ejected through interactions with the central binary stars. We find significant differences between planets ejected through planet-planet interactions and those by the binary stars. The main differences that arise are in the distributions of excess velocity, where binary stars eject planets with faster velocities. These differences should be observable amongst known free-floating planets in nearby star-forming regions. We predict that targeted observations of directly imaged free-floating planets in these regions should be able to determine their preferred formation pathway, either by planet formation in single or multiple stellar systems, or through processes akin to star formation. Additionally the mass distributions of free-floating planets can yield important insights into the underlying planet populations. We find that for planets more massive than 20 $\rm M_{\oplus}$, their frequencies are similar to those planets remaining bound and orbiting near the central binaries. This similarity allows for effective and informative comparisons between mass distributions from microlensing surveys, to those of transit and radial velocities. Ultimately, by observing the velocity dispersion and mass distribution of free-floating planets, it will be possible to effectively compare with predictions from planet formation models, and to further understand the formation and evolution of these exotic worlds.


INTRODUCTION
Until recently planets have only been observed orbiting either a single star (see Winn & Fabrycky 2015, for a review) or pairs of stars (e.g.Kepler-16b (Doyle et al. 2011) and recently BEBOP-1c (Standing et al. 2023)).In the last decade however, microlensing surveys have discovered a handful of planets that are not bound to any known stellar system (Mróz et al. 2019(Mróz et al. , 2020;;Koshimoto et al. 2023;Sumi et al. 2023).Such "free-floating planets" (FFPs) have now been found to exhibit a variety of planet masses, ranging from terrestrial (Mróz et al. 2020), to Neptunian (Mróz et al. 2019;Koshimoto et al. 2023) to Jovian (Sumi et al. 2011).Additionally, direct imaging of nearby clusters have also discovered numerous planetary mass objects that are not attached to close-by systems, whilst also showing evidence of planetary multiplicity (Pearson & McCaughrean 2023).With more of these planets now being discovered, their formation remains an interesting and unsolved question.
Whilst the mechanisms for the formation of FFPs are numerous, those involving planet formation have not explored the formation of FFPs from full population models, including the growth and evolution of the planetary systems.Additionally they generally focus on specific types of planet, e.g.giant planets, or by demonstrating efficient mechanisms for FFP formation.Previous global models of circumbinary planet formation showed that planets were frequently ejected from the systems, as they formed and grew in the circumbinary disc (Standing et al. 2023;Coleman et al. 2023Coleman et al. , 2024a)).These outcomes were all natural byproducts of the planet formation process where the models were able to form planetary systems similar to  In this paper, we utilise these models and run numerous simulations exploring the effects of different initial conditions on the planet formation process, and the final planetary systems.We focus on the formation and ejection of FFPs that arise naturally in the simulations, and find numerous trends within the distributions of planet properties.These include being able to differentiate between whether planets are ejected by the binary stars or through planet-planet interactions, and also the effects of different parameters on the mass, number and velocity distributions of FFPs.The distributions and predictions derived from our results will provide useful comparisons with present and future observations of FFPs, both those found through microlensing surveys finding old and evolved planets (Sumi et al. 2023), and additionally those found through direct imaging within nearby clusters (Pearson & McCaughrean 2023).
This paper is laid out as follows.We outline our physical model in Sect.2, whilst we describe our population parameters in Sect.3. In Sect.4, we outline the results of our population models, whilst in Sect. 5 we explore the effects that changing initial conditions have on the distributions of FFPs.Finally, we discuss our results and draw conclusions in Sect.6.

PHYSICAL MODEL
In the following section, we provide a basic overview of the physical model we adopt and the numerical scheme used to undertake the simulations.The N-body simulations presented here were performed using the mercury6 symplectic N-body integrator (Chambers 1999), updated to accurately model planetary orbits around a pair of binary stars (Chambers et al. 2002).We utilise the 'close-binary' algorithm described in Chambers et al. (2002) that calculates the temporal evolution of the positions and velocities of each body in the simulations with respect to the centre of mass of the binary stars, subject to gravitational perturbations from both stars and other large bodies.We also include prescriptions for the evolution of 1D protoplanetary discs as well as disc-planet interactions.With the disc models being 1D in nature, we also include prescriptions that take into account non-axisymmetric effects (i.e. a precessing eccentric inner disc cavity) due to the binary stars.Details on the full model and the additional prescriptions due to the binary can be found in Coleman et al. (2023Coleman et al. ( , 2024a)).However we briefly describe the model below.
(i) We solve the standard diffusion equation for a 1D viscous α-disc model (Shakura & Sunyaev 1973;Lynden-Bell & Pringle 1974).Disc temperatures are calculated by balancing black-body cooling against viscous heating and stellar irradiation from both stars.The viscous parameter α remains constant throughout most of the disc, but increases close to the central stars to mimic the eccentric cavity that is carved out by the tidal forces of the central stars.When giant planets are present, tidal torques from the planets are applied to the disc leading to the opening of gaps (Lin & Papaloizou 1986).
(ii) We incorporate models of photoevaporative winds removing material from the disc, both internally driven through radiation emanating from the central stars, and externally driven due to FUV radiation from nearby sources (i.e.O/Bstars).For internal photoevaporation, mainly due to EUV and X-ray radiation from the central stars (Clarke et al. 2001;Owen et al. 2010), we include the photoevaporation models of (Ercolano et al. 2021;Picogna et al. 2021) where the wind is assumed to be launched thermally from the disc upper and lower surfaces beyond a critical radius.We assume the most massive star dominates the high-energy radiation driving internal photoevaporation, and assume it has an X-ray Luminosity LX = 10 30.5 , consistent with the aveage value for Solar mass stars in star forming regions (Flaischlen et al. 2021).With ionising radiation not just impacting the disc from the central stars, but also from nearby stars in the local star-forming region (e.g.Haworth et al. 2018Haworth et al. , 2023)), we include external photoevaporation in the models to account for the effects of ionising FUV photons.We adopt the model found in (Qiao et al. 2023) which drives a wind from outside the radius where the disc becomes optically thin.
(iii) Eccentric cavities, arising because of the tidal torque from the central binary, have been seen in observations and numerical simulations of circumbinary discs (e.g.Artymowicz & Lubow 1994;Dutrey et al. 1994;Pierens & Nelson 2013;Mutter et al. 2017;Thun et al. 2017;Coleman et al. 2022).The shape and size of these cavities depends on the binary properties, i.e. mass ratio and eccentricity, and local disc properties such as the viscosity parameter α (Kley et al. 2019).The main effect of circumbinary disc cavities on planet formation is the creation of a planet migration trap as the corotation torque is increased due to the positive surface density gradient at the cavity edge.To simulate the effects of an eccentric cavity in our 1D disc model, we ran 2D hydrodynamical simulations of circumbinary discs using fargo3d (Benítez-Llambay & Masset 2016) to determine the cavity structure.In the 1D models, we simulate the azimuthally averaged surface density profile of the cavity by adjusting the viscosity parameter α, whilst maintaining a constant gas flow rate through disc.This forms an inner cavity in the disc and leads to a buildup of material at the outer edge of the cavity, as required.
(iv) Using the above 2D hydrodynamical simulations, we take into account the precessing, eccentric nature of the inner disc cavity in 1D models through construction of 2D maps of the gravitational acceleration experienced by test particles embedded in the disc due to the non-axisymmetric density distribution.We also create maps of the gas surface densities and gas velocities.These maps are used when integrating the equations of motion of planets and when calculating relative velocities between planets and drifting pebbles.These effects are mainly relevant near the cavity region.
(v) We include the planetesimal and embryo formation models found in Coleman (2021).We assume that as the pebble production front moves outward and pebbles drift inward, they can collect in short-lived traps due to non-axisymmetric perturbations in the discs.As solids collect in the traps, then planetesimals may form through gravitational collapse, when the local particle density exceeds the Roche density, and when the local dust-to-gas ratio exceeds unity (Johansen et al. 2007(Johansen et al. , 2009)).The size distribution of the planetesimals then follows a power law plus an exponential decay (Johansen et al. 2015;Schäfer et al. 2017;Abod et al. 2019).We assume that the most massive planetesimal that forms in each collapse event is ultimately a planetary embryo that is able to dominate the accretion and dynamical evolution of other planetesimals that formed in that vicinity.This planet is then able to accrete pebbles that are drifting past its orbit (Johansen & Lambrechts 2017), and planetesimals from the newly formed local reservoir (Fortier et al. 2013), whilst also undergoing mutual interactions and collisions with other planetary embryos.
(vi) The main source of the accretion of solids in our models is through pebble accretion.We follow the pebble accretion model of Johansen & Lambrechts (2017), where a pebble production front moves outwards in the disc over time.This production front arises from dust particles coagulating and settling to the disc midplane forming pebbles.Once these pebbles become large enough, they begin to drift inwards through gas drag forces, thus creating a pebble production front when the drift time-scale is equal to the growth time-scale.As the pebbles drift inwards they can be accreted by planetary embryos, allowing them to grow on short time-scales (Lambrechts & Johansen 2012).The accretion of pebbles continues until the planets reach the pebble isolation mass, that being the mass where planets are able to sufficiently perturb the local disc, forming a pressure bump exterior to the planet's orbit, that traps pebbles and halts pebble accretion on to the planet (Lambrechts & Johansen 2014;Ataiee et al. 2018;Bitsch et al. 2018).
(vii) The accretion of gaseous envelopes on to solid cores occurs once a planet's mass exceeds 1 M⊕.We utilise the formulae based in Poon et al. (2021) that are based on fits to gas accretion rates obtained using a 1D envelope structure model (Papaloizou & Terquem 1999;Papaloizou & Nelson 2005;Coleman et al. 2017).To calculate these fits Poon et al. (2021) performed numerous simulations, embedding planets with initial core masses between 2-15 M⊕ at orbital radii spanning 0.2-50 au, within gas discs of different masses.This allowed for the effects of varying local disc properties to be taken into account when calculating fits to gas accretion rates, a significant improvement on fits used in previous work (e.g.Coleman & Nelson 2014, 2016a,b).We use these fits until a planet is massive enough to undergo runaway gas accretion and open a gap in the disc.The gas accretion rate is then limited to either the maximum value of the fits from Poon et al. (2021), or the viscous supply rate.All gas that is accreted onto a planet is removed from the surrounding disc, such that the accretion scheme conserves mass.
(viii) We use the torque formulae from Paardekooper et al. (2010Paardekooper et al. ( , 2011) ) to simulate type I migration due to Lindblad and corotation torques acting on planetary embryos.Corotation torques arise from both entropy and vortensity gradients in the disc, and the possible saturation of these torques is included in the simulations.The influence of eccentricity and inclination on the migration torques, and of eccentricity and inclination damping are included (Fendyke & Nelson 2014;Cresswell & Nelson 2008).(ix) Type II migration of gap forming planets is simulated using the impulse approximation of Lin & Papaloizou (1986), where we use the gap opening criterion of Crida et al. (2006) to determine when to switch between type I and II migration.Thus, when a planet is in the gap opening regime, the planet exerts tidal torques on the disc to open a gap, and the disc back-reacts onto the planet to drive type II migration in a self-consistent manner.

POPULATION PARAMETERS
Previous work showed that ejected planets are a natural outcome of circumbinary planet formation simulations (Standing et al. 2023;Coleman et al. 2023Coleman et al. , 2024a)).Those studies varied the initial disc mass, metallicity, strength of the external environment, and the level of turbulence in the disc.However their main focus was on the types of circumbinary systems that form, and not on the properties of planets ejected from those systems.For the population of circumbinary systems in this work, we vary the initial disc mass, the binary separation, the strength of the external environment, and the level of turbulence in the disc.We only explore a single combined stellar mass, mass ratio and binary eccentricity, since for each combined set of those parameters, the properties of the central cavities are different, and so varying them would require large numbers of hydrodynamic simulations to be run to calculate the properties of the cavities.We base the central stars on that found in the BEBOP-1 circumbinary system (Kostov et al. 2020;Standing et al. 2023).The stellar and model parameters used to simulate the central cavities in the circumbinary discs as described above can be found in Table .1.The reason we choose the BEBOP-1 parameters for our simulations, is that the 1D circumbinary disc models require hydrodynamical simulations to be performed in order to provide prescriptions for the circumbinary cavity, disc eccentricity and the multi-dimensional maps for the resulting torques and gas velocities (see points iii and iv in Sect.2).These simulations were already performed in recent work (Coleman et al. 2024a).Additionally, choosing the BEBOP-1 parameters allowed preliminary comparisons with previous work that contained some overlapping simulations, but concentrated on different objectives, i.e. forming systems similar to BEBOP-1 (Coleman et al. 2024a).The choice of a total combined stellar mass of 1.3 M⊙ is also similar to the average combined stellar mass for binary stars (Raghavan et al. 2010, ∼ 1.5 M⊙).Table 2 shows the parameters that we do vary, and we take random values between the limits shown in table 2, with the last column denoting whether we randomise in log or linear space.
We choose initial disc masses between 5 and 15% of the combined binary mass, with the maximum disc mass being equal to the most massive disc a star can host before it becomes gravitationally unstable (Haworth et al. 2020).Numerous works have provided observational estimates for the viscosity parameter α (Isella et al. 2009;Andrews et al. 2010;Pinte et al. 2016;Flaherty et al. 2017;Trapman et al. 2020;Villenave et al. 2020Villenave et al. , 2022)), whilst theoretical work has contributed additional constraints (Standing et al. 2023;Coleman et al. 2024a).We adopt values of α that are consistent with such estimates (see Rosotti 2023, for a recent review).For the final parameter that affects the disc lifetime, the rate of external photoevaporation, we vary the strength of the local environment UV field, ranging from 1 G0 to 10 5 G01 .These values are consistent with what is expected across star forming regions such as Orion, with other low mass regions such as Taurus and Lupus occupying the lower region of our parameter space (Winter et al. 2018).The final parameter we vary is the binary separation, which we explore separations between 0.05 and 0.5 au.
We initialise the surface density in the circumbinary disc following Lynden-Bell & Pringle (1974) where Σ0 is the normalisation constant set by the total disc mass, (for a given RC), and RC is the scale radius, which sets the initial disc size, taken here to be equal to 50 au.The domain of our circumbinary discs extend from the binary separation as the inner edge to an outer edge of 500 au where for all discs the surface density is nearly always at our floor value of 10 −5 gcm −2 .We run each simulation for 10 Myr to account for the entire circumbinary disc phase, and additionally allowing for dynamical evolution of the systems after the dispersal of the circumbinary discs.Planets are removed from the simulations once they either collide with the central stars or other planets or when they are ejected from the system.We define a planet as having been ejected once it enters into a hyperbolic orbit, i.e. having an eccentricity greater than 1, and when its distance from the barycentric centre of the system has exceeded 1000 au.

RESULTS
The main objectives of this work are to explore the properties and distribution of free-floating planets that originate in circumbinary systems before being ejected.Whilst we do not discuss the properties of the remaining circumbinary systems themselves, or their formation as a whole, their formation pathways and outcomes are broadly similar to those found in Coleman et al. (2024a).We will explore the effects that the varying of the initial parameters in this work have on those outcomes in future work.Before we explore the distributions of ejected planets, and how the initial conditions affects their properties, we will describe the evolution of a typical simulation that led to multiple ejections from a single circumbinary system.

Example system with multiple ejection sites
The formation pathways of the planets in the circumbinary systems loosely follows that described in previous works (e.g.Coleman et al. 2023Coleman et al. , 2024a)).The main difference to those works here, is that by including the planetesimal and embryo formation models of Coleman (2021), the models form the initial planetary embryos, instead of starting with a collec-tion of embryos spread throughout the disc.These embryos were then able to start accreting pebbles and planetesimals from their local surroundings, allowing them to grow.Once the planets reached ∼lunar masses, then they were able to migrate in the disc, encountering other planets, with some mutual interactions leading to collisions or ejections.Additionally, as the planets reached terrestrial-super-Earth masses, they were able to accrete gaseous envelopes.The planets mainly migrated inwards, towards the central cavity, where the buildup of material at the cavity acted as a migration trap, allowing planets to congregate there, again inducing collisions and ejections.The architectures of the final planetary systems, and the planets found within them, were similar to those formed in Coleman et al. (2024a).
To highlight the evolution of a single simulated system, Fig. 1 shows the temporal evolution of planet masses (top panel), and semimajor axes (bottom panel).We highlight those planets that are ejected from the simulation and with masses mp > 1 M⊕ with coloured lines, whilst the grey lines show all other planets.The crosses show the final planet masses and semimajor axes from when before they were ejected, with black crosses showing planets ejected through interactions with the binary stars, and the red crosses being through planet-planet interactions.Looking at the left hand side of Fig. 1 it is clear that the planets are forming with masses mp = 10 −3 -10 −4 M⊕, over the first 0.1Myr with those closest to the star forming first, around ∼ 4 au, with the outermost planet forming at∼ 40 au.Over the next 0.2 Myr, the planets accrete pebbles and planetesimals, whilst undergoing mutual collisions with other planets.This allowed a number of planets to reach terrestrial and super-Earth masses, where they could then begin to migrate through the disc.This evolution is nicely highlighted by the lower green and dark red lines in Fig. 1 with the dark red line reaching a mass mp ∼ 20 M⊕ after 0.2Myr through accreting large amounts of pebbles and planetesimals and colliding with multiple growing embryos.
With multiple super-Earths and Neptune mass planets forming, this leads to a number of interactions with a few planets being ejected.As well as low mass planets being ejected, three planets with masses mp > 1 M⊕ are ejected around 0.3 Myr.This cluster of ejections is due to the interactions between the green and dark red planets orbiting closer to central binary around ∼ 1 au.As these planets grew, they underwent runaway gas accretion, becoming giant planets.They then interacted with each other, forcing the green planet to move on to a more eccentric orbit, and to be scattered out to ∼ 5-10 au.Once there, it interacted with other planets, forcing them to be ejected from the system.This is highlighted by those planets ejected at around 0.3 Myr.Over the next 0.2 Myr, that giant planet migrated back in towards the central binary, with the inner giant being ejected after interacting with one of the stars.Additionally, as the green planet migrated inwards, it interacted with the dark blue planet, forcing it to interact with the binary and be ejected from the system.As the giant planet represented by the green line migrated closer to the binary, the planet shown by the light red line underwent runaway gas accretion, causing another episode of dynamical instabilities and ejections from the system.This resulted in the cluster of ejections at ∼ 0.5Myr.With fewer planets now orbiting in the system, the dynamical interactions were quenched and so the evolution of the system then followed a more orderly pathway, with the final system containing a giant planet orbiting at 7.8 au, and a few sub-terrestrial and terrestrial mass planets orbiting between 22-117 au.
The ejections in the simulation described above were commonplace amongst the simulations presented in this work.Large frequencies of dynamical interactions, especially around the cavity close to the binary stars, acted to increase eccentricities and force planet to interact with more massive objects, e.g.giant planets or the binary stars themselves, and by ejected from the system.Additionally the example shows the abundance of ejections from circumbinary systems, whilst simultaneously showing the formation of the surviving circumbinary planets that have been seen in observations.

Population statistics
With Sect.4.1 outlining an example of the formation and ejection of multiple planets from a circumbinary system, we now discuss the properties and the distributions of the population of ejected planets as a whole.We begin by exploring the number of planets ejected over the course of each simulation.

Frequency of ejected planets
Recent work has estimated that the number of free-floating planets per star in the galaxy is 21 +23 −13 with masses mp > 0.33 M⊕ (Sumi et al. 2023).Previous works looking at a small range of parameters have also shown that between 3-9 planets are ejected from circumbinary systems around Kepler-16, Kepler-34 (Coleman et al. 2023), and BEBOP-1 (Standing et al. 2023;Coleman et al. 2024a).The numbers of planets ejected per system here is in broad agreement with that expected from observations, and from previous theoretical work.We find that each circumbinary system ejected 5 +2 −3 planets with masses greater than 1 M⊕ over the 10 Myr runtime of the system.
In Fig. 2 we show the cumulative distribution function for the number of planets ejected per system.The blue line shows the number of planets ejected per system for all planet masses, whilst the red and yellow lines show the number of planets ejected with masses mp > 1 M⊕ and 100 M⊕ respectively.As mentioned above the average number of planets ejected per system with masses mp > 1 M⊕ is around 5 planets.These averages change to 14 planets of all masses per system, and 0.6 planets of masses mp > 100 M⊕.In regards to the number of giant planets, this is consistent with the expectations from previous simulations around Kepler-16 and Kepler-34 (Coleman et al. 2023).When comparing to observations, these values are towards the lower end on their predictions.However the predictions from Sumi et al. (2023) is only based on few detections, hence the large uncertainty in their predictions, and so with more observed FFPs, by the Nancy Grace Roman Telescope (Spergel et al. 2015;Bennett et al. 2018) or the Rubin Observatory (Ivezić et al. 2019) for example, then the observed predictions for the number of FFPs per star would become more robust.
Interestingly whilst the number of ejected planets may be less than that predicted from observations, the total mass in planets ejected per system is larger than in observations.Based on FFPs found in microlensing surveys, Sumi et al. (2023) predict a total mass of 80 +73 −47 M⊕.From our simulations here, whilst we only eject 5 planets on average that have masses mp > 1 M⊕, the average total mass of planets ejected per system is equal to 147 +233 −84 M⊕.With 0.6 giant planets being ejected per system, they obviously contain the majority of the mass ejected per system.Their variance in masses also provides significant deviations of the total ejected mass.Nonetheless, the simulations show that more mass is ejected from the systems compared to observations.These differences could arise form the lack of observed FFPs yielding poorly constrained estimates, but additionally the simulations presented here were based on a central stellar mass of ∼ 1.31 M⊙, whereas those planets observed in microlensing surveys could have originated around a wide variety of stars.Running such a population of circumbinary planet formation simulations following a stellar IMF (e.g.Kroupa 2001) is not yet feasible computationally, but simulating such a population could yield estimates that are more consistent with future observations.

Temporal distribution of ejections
Recent observations of nearby star forming regions have found numerous planetary mass objects, down to a mass of ∼0.5 Jupiter masses, whilst their formation pathways are still an open question (Pearson & McCaughrean 2023).With planet and star formation still occurring in these regions, this makes the simulations presented here ideal to determine when planets are ejected.As discussed above, numerous planets are ejected over the 10 Myr simulation time, including the lifetime of the circumbinary discs.However, when the planets are ejected could be of use for observation surveys of nearby star forming regions to determine whether it would be expected to observe FFPs that originate in circumbinary systems.
To that end, Fig. 3 shows the cumulative distribution function of planets as a function of their ejection time.The blue line shows the distributions for all planets, whilst the red and yellow lines again show the distributions for planets mp > 1 M⊕ an mp > 100 M⊕ respectively.As can be seen from the red and yellow lines, it takes the systems at least ∼ 0.2 Myr to begin to see these planets ejected from the systems.This is due to the need to form the initial planetary embryos, allow them to accrete pebbles and gas, and then interact with other planets, or the central binary itself, leading to their ejection.This is not so much the case for lower mass planets that can be ejected much earlier in the disc lifetime, through interactions with the binary stars shortly after their formation.Once sufficient planetary growth has occurred, it is clear that the rate at which planets are ejected remains roughly constant for the next 1-2 Myr.This arises due to the different growth time-scales of planet in different regions of the disc, as well as the dynamical nature of the ejected planet's evolution that causes it to move into a situation where may interact with other objects and be ejected from the system.Therefore, Fig. 3 shows that the majority of the ejections from circumbinary systems occurs in the first few million years of their lives, and so it would be expected to directly observe more FFPs in younger star forming regions with ages 0.5-3 Myr (e.g.ONC), than in older regions (e.g.Lupus).Comparisons of the number of FFPs observed in these regions would provide valuable information on whether FFPs originate around circumbinary stars, or whether other processes, i.e. star formation processes, are responsible for these objects.If star formation was delayed, or prolonged over a long period of time, then this would extend the times at which FFPs would be observable in star forming regions.Therefore, observations of large numbers of FFPs in older clusters, e.g.Upper Sco, would provide evidence for prolonged or multiple epochs of star formation in the region.

Mass comparisons to bound planets
Whilst direct observations of FFPs in the Trapezium Cluster have only found planets down to a mass of 0.5 Jupiter masses (Pearson & McCaughrean 2023), microlensing surveys are now finding planets with super-Earth or even terrestrial masses (Mróz et al. 2019(Mróz et al. , 2020;;Koshimoto et al. 2023).With only a handful of sub-Jupiter mass FFPs found though, a robust estimation of the underlying mass distribution is not yet forthcoming, though future observations with the Roman Telescope or Rubin Observatory will bridge these uncertainties and deliver a mass distribution of FFPs.To that end, it is useful to understand the mass distribution of planets ejected from binary systems, and more importantly, how that distribution relates to the distribution of the remaining planets in the systems that can be detected through other detection methods (e.g.transit or radial velocity surveys).Understanding any differences in these populations will aid in the modelling of planet formation and evolution.
Figure 4 shows the cumulative distribution function of planet mass for those planets with masses mp > 1 M⊕.The blue line in Fig. 4 shows the distribution for planets ejected from the systems, whilst the red and yellow lines show the distributions for all bound planets, and for bound planets with semi-major axes ap ≤ 10 × a bin , i.e. with orbital periods amenable for detection by transit or radial velocity surveys.It is clear to see that for all three distributions the bulk of the planets that remain bound and ejected are of low mass with 70-80% of planets having masses mp ≤ 20 M⊕.There are differences there however, where there are few terrestrial planets ejected, and considerably more super-Earths and Neptune mass planets than for the distribution of the bound planets.For the bound planets the increase in the distributions are steady from terrestrial mass up to ∼ 15 M⊕ showing where pebble accretion dominates the formation of these planets.However, for the ejected planets, more super-Earth and Neptune mass planets are ejected since as multiple planets of these masses form in a single system, they can excite eccentricities and force other planets on to orbits where they interact with the binary stars and get ejected.Systems of Earth-mass planets are less likely to undergo this evolution since they need to be orbiting with significantly closer proximity in order to dynamically excite each other on to the required eccentric orbits for interactions with the binary stars.
Interestingly for planet masses greater than 20 M⊕, there is very little difference between the mass distributions of ejected planets (blue line) and those that remain bound to the star but orbit within 10 × a bin (yellow line).This shows that of the planets that form in those systems, similar fractions of planets as a function of planet mass remain bound and are ejected.The inset plot in fig. 4 shows the cumulative distribution function for planets with masses mp ≥ 20 M⊕ for ejected planets and those planets that remain bound but with semimajor axes ap ≤ 10 × a bin .It clearly shows the similarity between the two distributions for plants with masses greater than 20 M⊕.This result indicates that the observed mass distributions of FFPs through microlensing surveys should be similar to the mass distribution of planets seen in transit and radial velocity surveys.Obviously this would require a large number of observed planets, since the stellar populations that account for the observed planets would also need to be consistent with each other.Additional effects that affect the formation such as the turbulence in the disc, and the local star forming environment, which we will discuss later, can also affect the distributions of planets and so differences of the mass distributions in these planets could also indicate different initial conditions for such planets.An example of this could include comparing FFPs that formed and are found in the galactic bulge, to those planets that formed in nearby low-mass star forming regions.Ultimately though, the observed mass distributions of these Neptune to Jupiter mass planets can inform significantly on the formation history of such planets.

Velocity signatures
The mass distributions shown above are useful for comparing with FFPs found with microlensing surveys.However they are not as useful to compare with the more massive planets found with direct imaging by for example JWST or Roman.One property that is useful to compare with such observations is the observed velocity differences to other stars in the star forming regions.Numerous works have explored the velocity dispersions of stars in nearby star forming regions.For example the velocity dispersions found in the Orion Nebula Cluster is equal to σv ∼1-3 kms −1 (van Altena et al. 1988;Kim et al. 2019;Theissen et al. 2022).It would be expected that planetary mass objects that form at the tail end of the star formation process would have similar velocity dispersions.However if such objects formed through planet formation processes and were ejected from their systems by other planets or a binary star for example, then their velocity dispersions could be considerably different.
For an object to be ejected from a system, their velocity must be larger than the escape velocity.If such objects are launched on these hyperbolic orbits with more velocity than the escape velocity, then they will have an excess velocity which can be measured as these objects move through the local regions.To calculate the excess velocity we use the hyperbolic orbital elements of the planet as they are removed from the simulations when they reach 1000 au.More specifically it is equal to: where v, r, and a are the planet's velocity, distance from the barycentre, and semi-major axis when removed from the simulation.
In Fig. 5 we show the planet excess velocities as a function of the ejected planet mass.The blue points represent individual planets, whilst the black points show a binned average with a bin size equal to 0.05 dex.As can be seen in Fig. 5 there is a large variance in the excess velocity of objects, especially for those of lower mass.For lower mass planets, they can be deflected on to extremely close encounter orbits with one of the binary stars, which can result in extremely fast ejection velocities.This gives rise to the excess velocities for these objects reaching ∼100 kms −1 .Looking at the black points showing the average values, it can be seen that the average ejection velocities increases from 8 kms −1 for Earth mass planets up to 13 kms −1 for 15 M⊕ planets.This increase is due to multiple planets forming and mutually interacting close to the edge of the central cavity, allowing planets to strongly interact with the binary stars.Interestingly the variance in the excess velocity, and the average velocity decreases as the ejected planet mass increases.By looking at the blue points it is clear that the majority of planets with masses mp > 30 M⊕ are ejected from their systems with v∞ ≤20 kms −1 , whilst the average excess velocity falls from 12 kms −1 for 30 M⊕ planets down to 8 kms −1 for 500 M⊕ planets.When planets reach these masses, the damping forces they feel from the circumbinary discs act to reduce their eccentricities resulting in weaker interactions with the binary, whilst there are also fewer planets of a similar mass that can force them to interact with the binary through strong dynamical encounters.Nonetheless, it is interesting here that the average velocity of planets ejected form circumbinary systems is around 10 kms −1 , a factor few times larger than the velocity of other objects formed in the star forming regions.
Whilst Fig. 5 showed the excess velocity and masses of ejected planets, it is not able to answer the question if there are differences between those values for planets ejected through interactions with the binary stars or through interactions with planets further out in the system.The origins of the ejection event are important, since those that are ejected through interactions with other planets, far from the binary stars, will have properties akin to planets ejected in single stellar systems of a similar central mass.Dynamically, when not close to the binary, the escape velocity from mutual planetary interactions has to be larger than the escape velocity of the orbit, identical to the requirements in single stellar systems.Therefore it is possible to compare the properties of ejected planets from interactions with the binary, i.e. binary induced, to those where the binary is not involved akin to single stellar systems.In order to make this comparison we must first differentiate within our ejected sample those planets that are ejected through interactions with the binary, and those from further out in the disc.Using the hyperbolic orbital elements, we can calculate the pericentre of the orbit.Under the assumption that the last significant gravitational interaction was what led to the ejection of the planet, then the pericentre of the orbit is the approximate location of this event, thus informing us where the planet was ejected from.
Figure 6 shows the ejected planet mass versus the pericentre distance, i.e. the ejection location.The colours show the excess velocity of the planets.From Fig. 6 we can immediately identify three separate populations.The first population of planets are those with pericentre distances pp > 10 au, and masses, 10 −4 M⊕ < mp < 10 M⊕.These low mass planets and planetary embryos are ejected far from the binary stars, typically through mutual interactions with more massive objects, i.e.Neptune to Jupiter mass planets depending on their orbital separation.These are the planets that will be most like those that are ejected in single stellar systems, since the binary played a negligible role in the ejection of these objects.Interestingly the planets ejected in this region, typically have excess velocities around 1-10 kms −1 , showing that they are weakly ejected from the system, and would have similar velocity dispersions to surrounding objects in their local region.
The second population can be found closer to the central binary and encompasses the low mass planets with pericentre distances pp < 20a bin , and masses, 10 −5 M⊕ < mp < 10 M⊕.Previous works have shown that circumbinary discs are dominated through interactions with the central binary within 20a bin , resulting in the formation of a cavity and a trapping location for migrating planetesimals, planetary embryos, and planets (Coleman et al. 2022(Coleman et al. , 2023)).This allows objects to concentrate in this region before being excited on to orbits that interact with the central binary, and thus then ejected.The effect of the binary is also extremely clear here on the excess velocity, with those objects that are ejected closer to the cavity edge, located typically around 1 au being ejected at much larger velocities, shown by the orange and yellow points in Fig. 6.Such ejections would not be possible in single stellar systems since the objects are too deep in the gravitational wells to be ejected by most planets.
The third and final population is those planets with masses mp > 10 M⊕.These planets have a wide range of pericentre distances from 0.01 au-1000 au.They are also typically ejected with excess velocities between 3-30 kms −1 .These ejection locations and velocities are very different to the two other populations, especially the pericentre distances.This arises due to numerous complications in calculating the ejection location.Firstly, the damping that acts on these planets can result in subtle changes to the planet's velocity and thus affect the excess velocity derived from the planet's orbital elements when it reaches 1000 au.The other main complication arises through the fact that these simulations involve multiple planets, and typically when planets are excited on to eccentric orbits, they are done so by sufficiently massive bodies.As planets are placed on to hyperbolic, escaping, orbits by the binary stars, there is a non-negligible possibility that these planets interact with other massive objects as they exit the inner parts of the system.Such interactions can again act to speed up or slow down the ejecting planet, sometimes causing it to be bound to the system again.Nevertheless these interactions add significant contamination to determining where the planets are ejected from.Aside from those complications, it is still important to acknowledge that these more massive planets have significantly large excess velocities, much greater than what is see for stellar populations in star forming regions.Such large differences in velocity dispersions of FFPs if shown would indicate that they originated in circumbinary systems.
To further the question of determining if it is possible to tell if FFPs originate in circumbinary systems, we plot the cumulative distribution function of the ejected planets excess velocities in Fig. 7.We only include planets with masses mp > 1 M⊕ to be more consistent with what may be observable with future observations.We show the CDF for all planets with the yellow line, whilst we split the population into two groups, those that are scattered through planet-planet interactions similar to the first population described above but including more massive planets (blue line), and those that are ejected through interactions with the binary stars (red line).For the planets that have complicated histories, i.e. those that may have interacted with other objects on their escape trajectory, we check whether their semi-major axis within the last 10,000 years of their time bound in the simulations was situated within 10a bin of the central binary.It is clear from Fig. 7 that there is a large difference in the velocity dispersions between those planets ejected through planet-planet interactions and those through interactions with the binary itself.Indeed, for planet-planet interactions, the average velocity that planets were ejected at was equal to 4.54 +1.21  −2.58 kms −1 , with the limits showing the interquartile range.This is in stark contrast to those ejected through interactions with the binary where the average velocity was equal to 11.97 +3.37  −5.54 kms −1 .Figure 7 therefore shows that the planets that are ejected through interactions with the binary stars have velocities roughly three times larger than those through planet-planet interactions, a result which should also hold when comparing ejections from binary systems to single stellar systems.With such large velocities, it would also be expected that these planets would have insufficient time to standardise their ve- locities to that of the stars in the star forming region.This is due to the ejection times of planets in most star forming regions (e.g.0.5 Myr for a region of width ∼few parsec) being much shorter than the relaxation time (Wang et al. 2015).Additionally, with such a large velocity signature, compared to the background stars in the local region, finding such a large velocity dispersion would be a clear signature of circumbinary planet formation, and ejection of FFPs.
The final point to take from Fig. 7 is the frequency at which ejections occur for both populations.As can be seen there is very little difference in the velocity dispersions between the binary and combined lines in Fig. 7.This shows that the binary induced ejections are contributing the most to the observed velocity dispersions.Indeed, when looking at the absolute numbers of ejections, binary induced ejections contribute 99.4% of all ejections from the system, highlighting how difficult and rare it is for planet-planet interactions to lead to ejections of an object.From observational surveys, 5-10% of Solar type stars in binary systems are of configurations comparable to those explored here, i.e. close-binaries (Offner et al. 2023).Combining this with Solar-type stars being found to contain 0.6 companions per stars (Offner et al. 2023), then for every close binary system explored here, there would be 7-14 singular Solar type stars.Coupling this with the numbers of ejected planets per system described previously, and the predicted numbers of FFPs per star based on microlensing surveys, this would again indicate that binary systems, instead of single stellar systems, are responsible for the formation of FFPs.
Additionally, stellar flybys could be an efficient production mechanism of FFPs (Wang et al. 2023).However, the frequency of stellar flybys is determined by the number density of stars in a star forming region.For nearby star forming regions, e.g.ONC, the number density reaches 10 4 pc −3 (Winter et al. 2018), which equates to a background flux of 10 4 G0 in this work.For this number density, the closest encounter for stellar flybys lies between 100-1000 au, much too far from the inner system to drive efficient ejections of planets (Winter et al. 2018).For less dense clusters, the closest approaches are Confirming some of the other observational signatures, such as the velocity dispersion would further add confidence to the conclusion that binary stars are essentially dominant FFP factories.

HOW DO SIMULATION PARAMETERS AFFECT THE DISTRIBUTIONS
Section 4 explored the population statistics as a whole that arose from the simulations.We now determine the effects that different initial parameters have on the observable distributions of FFPs.These include the level of turbulence in the disc and the external photoevaporation rate, which have previously been found to significantly affect the evolution of protoplanetary discs and the planets that form within them (Coleman & Haworth 2022;Coleman et al. 2023).We will also explore the effects of the binary separation on the distributions.

How turbulence affects the FFP population
The level of turbulence in the disc is determined by the strength of α where in our simulations we have used values between 10 −4 -10 −2.5 , consistent with observations (Rosotti 2023).In Fig. 8 we show three separate plots similar to Fig. 6 where each plot shows the planet mass versus the pericentre distance (ejection location), with the colours showing the excess velocity of the ejected planets.The left-hand panel shows for planets that formed in discs with 10 −4 ≤ α < 10 −3.5 , the middle panel being 10 −3.5 ≤ α < 10 −3 and the right-hand panel being for planets in discs with 10 −3 ≤ α < 10 −2.5 .The different panels essentially show the effects of increasing turbulence from left to right.
Comparing the panels in Fig. 8 it is clear that there are differences depending on the level of turbulence.The main differences arise in the outer regions of the disc where there Figure 9. Cumulative distribution functions for number of planets ejected from each system.The level of turbulence in the disc is denoted by the colours with the blue line showing discs with 10 −4 ≤ α < 10 −3.5 , the red line showing discs with 10 −3.5 ≤ α < 10 −3 , and the yellow line being for discs with 10 −3 ≤ α < 10 −2.5 .are fewer planets ejected through planet-planet interactions as α increases (including sub-terrestrial mass planets).Additionally, fewer planets are ejected through interactions with the central binary.The main causes for these reductions is on the effectiveness of forming the initial planetary embryos and planetesimals decreasing as α increases.As α increases, the turbulence in the disc stirs pebbles to higher altitudes, resulting in a less dense pebble layer in the midplane.Since one of the criterion for the gravitational collapse of a pebble cloud to occur is that the local solids-to-gas ratio must exceed unity, the reduction in density due to higher α values makes this harder to achieve.Thus fewer planetesimals and planetary embryos are able to form and grow, and so there are then fewer planets interacting with each other or the binary and ejected from the system.
Figure 9 shows the number of all planets that are ejected from discs with different levels of turbulence, with blue lines showing low α values, and red and yellow lines showing intermediate and high values respectively.The reduction in the number of planets ejected as α increases is clear to see, with discs evolving with low α values ejecting ∼ 3 times more planets on average than those with large α values.These trends show that if observations can accurately constrain the number of FFPs per star, then they can also provide insights into the fundamental properties of protoplanetary discs.
As figs.8 and 9 showed that fewer planets were ejected from systems as the viscosity parameter α increases, this also has an effect on the mass distributions of planets that are ejected.Figure 10 shows the mass distributions for planets above an Earth mass ejected from discs with different strengths of α, as shown by the different coloured lines.The solid lines show those planets ejected through interactions with the central binary, whilst the dashed lines show those planets lost through planet-planet interactions.An interesting feature that arises from Fig. 10, especially for the binary induced ejections, is that as the strength of turbulence α increases, the mass distribution of ejected planets shifts to more massive planets.Indeed in discs with strong levels of turbulence, α ≥ 10 −3 , 46 % of planets ejected have masses m > 100 M⊕, whereas in discs with weak levels of turbulence, giant planets make up 4% of those ejected.Therefore more lower mass planets are ejected in weaker turbulent discs.The main cause of this change in distributions as a function of α is due to planets being able to more easily form in low turbulent discs, as the pebbles are more confined to midplane of the disc, aiding both planetesimal/embryo formation, as well as pebble accretion rates.As planets reach the pebble isolation mass, their accretion rate slows since gas accretion is inefficient there.This results in numerous terrestrial and super-Earth mass planets occupying similar areas of space, mutually interacting and driving up eccentricities.They then interact with the binary and are ejected from the system.In discs with higher α, this is not so much the case, since there are fewer planets that initially form and grow, and so a greater percentage of planets can undergo runaway gas accretion, becoming giant planets, before they are ejected from the systems.
Interestingly, the distributions in Fig. 10 appears to become more bi-modal as α increases, with super-Earth-Neptune mass and giant planets making up the majority of planets ejected in high α discs.The dearth of planets with masses between ∼ 30 M⊕ and ∼ 100 M⊕ is a signature of runaway gas accretion, since one the mass of the gaseous envelope becomes comparable to the core mass, the self gravity of the envelope allows it to quickly contract, driving up accretion rates (Coleman et al. 2017).This runaway gas accretion halts when the supply of gas to the envelope becomes limited since the planet opens a gap in the disc.However in low α discs, since the majority of planets ejected are only of super-Earth to Neptune mass, since there is greater frequencies of dynamical interactions, the signature of runaway gas accretion is effectively diminished.Additionally, with low values of α, planets are more easily able to open gaps in the discs.Observing such a bimodality in observations would therefore then give hints as to the underlying properties of protoplanetary discs.

What is the role of the local stellar environment
Previous studies have shown that the local environment influences how planets form (Winter et al. 2022;Qiao et al. 2023) and therefore also on the resulting planet populations (Standing et al. 2023;Coleman et al. 2024a).Those works found that as the strength of the local environment increased, the possibility for more massive planets to form and grow diminished.This was due to the reduction in accretable solids since external photoevaporation from the nearby massive stars effectively truncated the disc to small sizes.Within our simulations presented here, we find similar effects on the planets that form, and thus on the final planetary systems.There are also noticeable effects on the properties of free floating planets that escape from systems in different environments, of which we will now discuss.
Similar to Fig. 8, Fig. 11 shows three separate plots of ejected planet mass versus the pericentre distance (ejection location) with the colours showing the remaining excess velocity.Instead of the different panels showing the strength of turbulence in the disc, this time they show planets that formed in different environments.Again, this increases from left-to-right with the left-hand panel showing planets that formed in weak UV environments (< 10 2 G0), the middle panel showing for intermediate environments (10 2 G0-10 4 G0), and the right-hand panel being for planets forming in strong UV environments (> 10 4 G0).When moving form left to right on the plot, it is clear that the left-hand and middle panels are similar in terms of the populations of planets at different ejection locations and of specific planet masses, but there are significant differences for the right-hand panel, the strongest UV environments.Most noticeably is the reduction in giant planets ejected at the top of the panel, with the majority of giant planets this time being ejected from around the edges of the circumbinary cavities.Both, the reduction of giant planets ejected, and the confinement of the ejection location, is a result of the reduction in solids available for accretion.With less solids, fewer giant planets were able to form, allowing more giant planets to remain in stable orbits and having fewer interactions with the binaries that led to their ejections.This is also the case for less massive planets, e.g. of super-Earths and Neptune mass.With fewer massive planets forming, there are fewer opportunities for the planets to be forced on to eccentric orbits before being ejected from the system.With Fig. 11 showing fewer planets being ejected from discs evolving in strong UV environments, we now look at the total mass of planets ejected from systems in different UV environments.Figure 12 shows the cumulative distributions of the total combined mass of planets ejected in weak (< 10 2 G0, blue line), intermediate ((10 2 G0-10 4 G0), red line), and strong UV environments (> 10 4 G0, yellow line).Further highlighting the similarities between the left-hand and middle panels in Fig. 11, the red and blue lines in Fig. 11 are effectively identical.This shows that there is minimal effect on the planets forming and being ejected in weak and intermediate UV environments (< 10 4 G0).Whilst this may show that weak and intermediate environments do not affect the formation of planets, it is worth noting that the protoplanetary discs they form in are affected by such weak environments, since the discs are truncated down to ∼few hundred au (Haworth et al. 2023;Coleman et al. 2024b).Whilst the weak and intermediate environments showed few differences in the total mass of planets ejected, those planets forming in stronger environments are on average less massive.Looking at the yellow line in Fig. 12, 36% of systems in strong UV environments ejected planets with a combined mass of > 100 M⊕.Comparing this value to the red and blue lines, this percentage increases to 57% and 58% respectively, further highlighting the increase in total planet mass ejected from those systems.Going to much less massive planets, 14% of systems in strong UV environments ejected planets with a total mass below 1 M⊕, compared to 5% in weak and intermediate environments.
As Figs. 11 and 12 show there are fewer planets ejected in strong UV environments, this should be an observable difference in observed star forming regions.With Fig. 3 showing most ejections occur between 0.5-3 Myr, observations of young star forming regions with a variety of different UV field strengths, then there should be more planets, and larger total combined masses in weaker UV regions.Interestingly there is a negligible difference in the mass distribution of the planets ejected themselves, but there are just fewer planets that form and are ejected.

Possibility of JUMBOs?
Recent observations of the ONC cluster by JWST found numerous planetary mass objects, appearing unbound from nearby stellar systems (Pearson & McCaughrean 2023).Interestingly, a significant fraction (∼ 10%) of these free floating objects were found to be be of binary nature, and were termed as "Jupiter-mass Binary Objects" (JUMBOs).These objects with masses 0.5MJup ≤ mp ≤ 14MJup, are typically found with separations asep < 300 au.Whilst in our simulations, 18% of systems ejected at least 2 planets with masses mp > 100 M⊕, it is extremely difficult to naturally form giant planets as binary objects, or even on coorbital orbits.Additionally it is even more unlikely to be able to eject both planets at a similar time, allowing them to remain gravitationally bound to each other.
By comparing the ejection time of planets in the same system, and additionally the distance between them, we can determine when planets are ejected either as binary objects, or fortuitously along similar trajectories.Unfortunately, nearly all of the giant planets ejected in our simulations leave the systems along vastly different trajectories.Indeed, the closest "pair" of giant planets ejected in the simulations, were separated by ∼ 990 au when the second planet was ejected from the simulation and was 1000 au from the central binary.Further showing the point of the large distances between planets ejected from the same system, looking at "pairs" of planets with masses mp > 10 M⊕, only 0.6% of them had separations less than 1000 au.This shows how extremely rare it is for binary objects to be ejected from within circumbinary systems, pointing to other mechanisms for the formation of the observed JUMBO population.Such other mechanisms could include, formation in stellar systems and fortuitous ejection through flybys of nearby stars (Wang et al. 2023) or through star formation pathways (Portegies Zwart & Hochart 2023).

Effects of changing binary separations
The final key input parameter that we test to explore the differences in the distributions of FFPs is the binary separation.Typically most works exploring the formation of circumbinary planets only use a specific set of binary parameters, including the binary separation (e.g.Coleman et al. 2023Coleman et al. , 2024a)), since they are usually examining the formation of specific systems.Whilst we do not vary the binary mass ratio or eccentricity, since they affect the prescriptions for the central circumbinary cavity (Thun et al. 2017;Kley et al. 2019), we do vary the binary separation between 0.05-0.5 au.Varying the binary separation has previously been unexplored in the formation of circumbinary planets, of which we will explore in future work.Here however, we examine the effects that the binary separation has on the properties and distributions of FFPs.
In Fig. 13 we again show three plots of planet mass versus the pericentre distance (ejection location) for planets forming in systems with different binary separations.The colours again show the excess velocities as the planets are ejected from the systems.The left-hand panel shows planets that formed in systems with separations 0.05 au ≤ a bin < 0.2 au, whilst the middle shows for 0.2 au ≤ a bin < 0.35 au, and the right-hand panel for 0.35 au ≤ a bin ≤ 0.5 au.Comparing the different panels, they all generally contain similar populations of planets, however there are some small and subtle differences.For the smallest binary separation (left-hand panel), there is a larger number of planets with masses between 1-100 M⊕ ejected from around the cavity edges, situated between 0.1-1 au.Comparing this to the right-hand panel for the most separated binaries, the area of ejection is much more confined to a smaller region.Additionally, when comparing the excess velocity of these planets, it is clear that those ejected from the closer binaries in the left-hand panel, have larger velocities than their wider binary counterparts.This is not unsurprising since the orbital and escape velocities involved there are much larger, and so any deviations on to eccentric orbits, where the planets more freely interact with the stars, can lead to more energetic escapes.Furthermore, with the binaries being more compact, and ejections occurring in closer proximity, the likelihood of a planet on a hyperbolic orbit interacting with both binary stars, as well as other planets in the vicinity, also increases, which gives rise to a larger spread in the pericentre distance, since their velocities will be altered from the original ejecting event.
Whilst the distributions of planet mass, number of planets ejected, and the total mass ejected, all remain similar as a function of the binary separation, there is one distribution that is significantly affected.As noticed in Fig. 13 the excess velocities that planets retained were larger for closer binaries than for wider binaries.Similar to Fig. 7, we plot the cumulative distributions of excess velocities for ejected planets with masses mp > 1 M⊕ in Fig. 14.Here we separate those planets ejected from the binary with solid lines to those through planet-planet interactions with dashed lines, whilst the colours show planets ejected from: close binaries (blue), intermediate binaries (red) and wider binaries (yellow).Whilst there appears to be little difference for those planets ejected through planet-planet interactions (though there are few planets ejected in this manner), there is a noticeable trend for those planets ejected through interactions with the binary stars.As the separation of the binary stars increases, the excess velocity distributions moves to lower values, approaching the planet-planet induced velocities.This would suggest that for even wider binaries, i.e. with separations of ∼few au, there would be little difference in the excess velocity distributions irrespective of the ejection mechanism.This is unsurprising, since the escape velocity for the wider orbit planets would be reduced, and additionally, strong enough interactions with one of the binaries to eject the planet would be more easily attained, without the planet having extreme close encounters with one of the stars.With close binaries being observed frequently near to the sun (Raghavan et al. 2010) and with observing campaigns finding large numbers of eclipsing binaries (e.g.Kepler, Prša et al. 2011), signatures in excess velocities should be readily observable in FFPs found in star-forming regions.

DISCUSSION AND CONCLUSIONS
In this work we have explored the formation of free floating planets (FFPs) within circumbinary systems, where the planets are initially able to form similarly to remaining bound planets, before they are ejected either through interactions with other planets or with the central binary stars themselves.We used an updated version of the N -body code mercury6 including the effects of a central binary, and coupled to this a self-consistent 1D viscously evolving disc model containing prescriptions for planet migration, accretion of gaseous en-velopes, pebble accretion, disc removal through photoevaporative winds and prescriptions taking into account the effects of the central binary such as an eccentric cavity.Within the resulting populations we derived distributions and make observable predictions that can be tested with future observations.We explored the effects of numerous initial parameters on the distributions of FFPs including, the strength of the local radiation environment, the level of turbulence in the disc, and binary separation.The main results from our study can be summarised as follows.
(1) In agreement with previous works (e.g.Coleman et al. 2023Coleman et al. , 2024a)), circumbinary systems are efficient factories of FFPs.On average, each binary system ejects between 2-7 planets with masses mp > 1 M⊕.When considering giant planets only, i.e. those with masses mp > 100 M⊕, then only 0.6 planets are ejected per system.Whilst the number of FFPs formed is consistent with other planet formation studies, they are reduced compared to preliminary observations in mircolensing surveys (Sumi et al. 2023).
(2) Most of the FFPs that form within circumbinary system are ejected between 0.4-4 Myr, typically whist there are still circumbinary discs.This is important when attempting to observe FFPs in young star-forming regions, since they should be abundantly found there.
(3) Comparing the mass distributions of ejected planets to those remaining bound and orbiting close to the central stars, we find that for planets with masses mp > 20 M⊕, the distributions are similar.This means that observed mass distributions and occurrence rates from microlensing surveys should be comparable to the true populations of planets in circumbinary systems, close to the stars, and additionally to other observation surveys, e.g.transit or radial velocity, of circumbinary systems which are generally biased towards planets orbiting close to the central stars.
(4) Along with the mass distributions of FFPs, another main observable has arisen within the population of FFPs.As the planets are ejected from the systems, they retain significant excess velocities, between 8-16 kms −1 .This is much larger than observed velocity dispersions of stars in local star forming regions (van Altena et al. 1988;Kim et al. 2019;Theissen et al. 2022), and so determinations of the velocity dispersions of FFPs in stellar clusters should differentiate whether those planets arose from either star or planet formation process.
(5) Additionally the velocity dispersions of FFPs ejected through interactions with the binary stars is ∼3 times larger than those ejected through planet-planet interactions.This means that observed velocity dispersions can determine whether binary or single star systems are the dominant origin of FFPs.From our simulations, we find that planet-planet interactions only account for 0.6% of all FFPs larger than 1 M⊕, indicating that binary systems are the most likely origin of FFPs when taking observational constraints on binary fractions into account (Offner et al. 2023).
(6) The initial parameters studied also affect the distributions of the properties of FFPs.Notably, the level of turbulence affects the number of planets ejected, and their accompanying mass distributions.Weaker level of turbulence favours larger numbers of planets ejected, with more less massive planets (96% are less massive than 100 M⊕), with stronger turbulence favouring more massive planets (46% with masses greater than 100 M⊕).The strength of the lo-cal environment generally affects the total mass of planets ejected from circumbinary systems, with those in strong UV environments being most severely affected.Both of these effects of turbulence and the environment, are due to the reduction in solid material available for accretion by planets, either through faster dispersal of the disc, or through extended pebble and dust layers that hinder planetesimal formation and pebble accretion.
(7) Finally, the separation of the binary is found to have little effect on the frequency of FFPs, and their resultant mass distributions.This is due to the main planet formation processes occurring far from the influence of the binary stars, and then migration transporting the planets to the vicinity of the binaries.The binary separation does however, affect the velocity distributions of ejected planets, with wider binaries causing a reduction in the average excess velocity since planets are more easily ejected without having to be deflected onto orbits that pass close to the individual central stars.
The simulations here show that it is possible to use observations of the distributions of FFPs to determine their origins, as well as further properties of the environment, both the discs and the star forming regions they formed in.Differences in the distributions of FFP masses, their frequencies, and excess velocities, can all indicate whether single stars or circumbinary systems are the fundamental birthplace of FFPs.However, whilst this work contains numerous simulations, and explores a broad parameter space, it does not constitute a full population of forming circumbinary systems.Such a population would include varying combined stellar masses, mass ratios, and binary eccentricities, that would be motivated by observations.Since the modelling of central cavities in circumbinary discs require specific prescriptions, depending on the binary properties, it is not currently feasible to derive such a population.Should such a population be performed in future work, then comparisons between that population and observed populations would give even more valuable insight into the formation of these intriguing objects.
In addition, planet formation models are constantly evolving and adding new physics.With the population of FFPs originating in different regions of the circumbinary discs, then they would accrete material with varied compositions, both in the gas, and in the solids.In future work we will incorporate compositional models (e.g.Thiabaud et al. 2014Thiabaud et al. , 2015) ) that will allow us to determine if there are any chemical signatures within the FFP populations that can point to not only the stellar environment in which they formed, but more tightly constrain the regions within their natal protoplanetary discs.Should spectroscopic observations be performed for known directly imaged FFPs, and should we know their origins, then this would greatly inform on planet formation models as a whole, and indicate where improvements are needed within our models.Only then, will we be able to develop a full and complete model of planet formation that can explain the populations observed today, both the FFP population, and those exoplanets found through other methods, e.g.transit or radial velocity surveys.

Figure 1 .
Figure 1.Temporal evolution of planet semimajor axes (bottom panel), and masses (top panel) for an example system that ejected multiple planets as described in sect.4.1.The coloured lines show those planets that are ejected from the system and with final masses mp > 1 M ⊕ .The grey lines show the tracks for all other planets.The black and red crosses show the final semimajor axes and masses of planets ejected through interactions with the central binary, and through planet-planet interactions respectively.

Figure 2 .
Figure 2. Cumulative distribution functions for the number of planets ejected from the simulated circumbinary systems.The distributions are shown for all planets (blue line), for planets with masses mp > 1 M ⊕ (red line), and for giant planets with masses mp > 100 M ⊕ (yellow line).

Figure 3 .
Figure 3. Cumulative distribution functions for the time that planets are ejected from circumbinary systems.The distributions are shown for all planets (blue line), for planets with masses mp > 1 M ⊕ (red line), and for giant planets with masses mp > 100 M ⊕ (yellow line).

Figure 4 .
Figure 4. Cumulative distribution functions for the mass of surviving planets with masses mp > 1 M ⊕ .The distributions are shown for all planets ejected from the systems (blue line), for all remaining bound planets (red line), and for all bound planets with semi-major axes ap ≤ 10a bin .The inset plot shows the cumulative distribution function for planets with masses mp > 1 M ⊕ , highlighting the similarity between ejected planets and bound planets with semi-major axes ap ≤ 10a bin .

Figure 5 .
Figure 5. Excess velocities of ejected planets as a function of planet mass.Blue points show individual planets, whilst black points show the average excess velocity within mass bins covering 0.05 dex.

Figure 6 .
Figure6.Planet mass versus the pericentre distance for ejected planets.The pericentre distance gives an approximate location for the final interaction that led to the ejection of the planet.The colour coding shows the excess velocities of the planets after they have been ejected.

Figure 7 .
Figure 7. Cumulative distribution functions (bottom panel) and probability distribution functions (top panel) of planet excess velocities for ejected planets with masses mp > 1 M ⊕ .The lines differentiate between planets ejected through planet-planet interactions (blue line), planets ejected through interactions with the binary stars (red line), and a combined distribution (yellow line).

Figure 13 .
Figure 13.Same as Fig. 8, but the panels show for planets ejected in systems with small binary separations (0.05 < a bin < 0.2 au, left-hand panel), intermediate binary separations (0.2 < a bin < 0.35 au, middle panel) and large binary separations (0.35 < a bin < 0.5 au, right-hand panel).

Figure 14 .
Figure14.Cumulative distribution functions of planet excess velocities for ejected planets with masses mp > 1 M ⊕ .The lines differentiate between planets ejected through interactions with the central binary (solid lines), and through planet-planet interactions (dashed lines).The colours show for planets ejected from systems with small binary separations (0.05 < a bin < 0.2 au, blue lines), intermediate binary separations (0.2 < a bin < 0.35 au, red lines) and large binary separations (0.35 < a bin < 0.5 au, yellow lines).

Table 2 .
Values for the parameters varied amongst the populations.