EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 688 (August 2024) A&A, 688 (2024) A128 Full HTML
Open Access
Issue
A&A
Volume 688, August 2024
Article Number A128
Number of page(s) 18
Section Astrophysical processes
DOI https://doi.org/10.1051/0004-6361/202449421
Published online 12 August 2024

© The Authors 2024

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model.

Open Access funding provided by Max Planck Society.

1. Introduction

Circumbinary disks are pervasive in the Universe. They exist across astrophysical scales, ranging from supermassive black hole (SMBH) binaries (Begelman et al. 1980; Haiman et al. 2009; Roedig et al. 2011; Franchini et al. 2021; D’Orazio & Charisi 2023) to stellar systems. For stellar binaries, such disks are present in the early stages during star formation (Dutrey et al. 1994; Mathieu et al. 1997; Takakuwa et al. 2014; Lacour et al. 2016; Tofflemire et al. 2017; Long et al. 2021) as well as during the late stages for post-asymptotic giant branch (AGB) stars (Deroo et al. 2007; van Winckel et al. 2009; Bujarrabal et al. 2013; Gallardo Cava et al. 2021) and in post-common envelope systems (Kashi & Soker 2011; Reichardt et al. 2019; Röpke & De Marco 2023; Tuna & Metzger 2023). They may further occur in triple systems when the outer star donates mass to the inner system (Di Stefano 2020; Hamers et al. 2022; Dorozsmai et al. 2024). One of the primary goals in the study of circumbinary disks is to understand how the interaction with the disk influences the orbital evolution of the binary. The intricate interplay between the time-varying binary potential, the disk’s morphological response, and the resulting forces on the binary orbit can lead to a wide array of outcomes including secular changes in the binary orbital elements, a variation in the mass accretion, and possibly even mergers.

A substantial body of research, spanning from early analytical and numerical investigations (Artymowicz 1983; Artymowicz et al. 1991; Lubow & Artymowicz 1996; Pringle 1991) to more recent hydrodynamic simulations (e.g., MacFadyen & Milosavljević 2008; Cuadra et al. 2009; Shi et al. 2012; D’Orazio et al. 2013; Muñoz et al. 2019), reveals the sensitivity of these outcomes to the parameters of the system. The current understanding delineates a complicated landscape within a parameter space that remained, until recently, relatively uncharted. The evolution of the orbital parameters is a nontrivial function of mass ratio (Ragusa et al. 2020; Duffell et al. 2020; Siwek et al. 2023a,b; Dittmann & Ryan 2024) and eccentricity (Miranda et al. 2017; Zrake et al. 2021; Siwek et al. 2023b); and can also be influenced by the disk aspect ratio (Young & Clarke 2015; Ragusa et al. 2016; Tiede et al. 2020; Dittmann & Ryan 2022, 2024), viscosity (Duffell et al. 2020; Dittmann & Ryan 2022, 2024; Turpin & Nelson 2024), disk-orbit inclination (Moody et al. 2019; Smallwood et al. 2022; Tiede & D’Orazio 2024; Martin et al. 2024; Bourne et al. 2023), and cooling rate (Wang et al. 2023a,b).

Because of these complex dependencies, fully characterizing binary-disk interactions requires extensive parameter exploration. These studies are typically performed either in 3D smoothed particle hydrodynamics (SPH, Lucy 1977; Gingold & Monaghan 1977) where the orbit of the binary is integrated self-consistently from the forces of the gas particles, or with grid-based simulations where it is typical to fix a Keplerian orbit and compute the change in orbital elements in post-processing from the gas forces. These studies can be computationally demanding because of the high spatial resolutions and long integration times required to achieve converged results, limiting the possibility of self-consistently following the long-term evolution of the system.

Due to the resources required by hydrodynamic simulations, previous works interested in modeling the long-term effect of the disk on the orbit preferred to rely on analytical prescriptions. In some cases, the employed prescriptions are based on analytical disk models, which did not take into account the complex structure and dynamics of the gas inside the central cavity (Rafikov 2016; Siwek et al. 2020; Izzard & Jermyn 2023; Wei et al. 2024) and sometimes neglected accretion forces (Dermine et al. 2013; Vos et al. 2015). Other works chose to arbitrarily fix the speed of orbital shrinking (Tokovinin & Moe 2020; Tuna & Metzger 2023) or neglected the angular momentum exchanged with the disk (Di Stefano 2020; Hamers et al. 2022; Dorozsmai et al. 2024).

In this work, we propose a new approach to compute the long-term interaction of the binary with the circumbinary disk. Our approach is based on the results of high-resolution hydrodynamic simulations, which take into account both the gravitational and accretion forces. Specifically, we integrated across the parameter space the derivatives obtained from hydrodynamic simulations in Siwek et al. (2023a,b), which we then compared to the results of Zrake et al. (2021) and D’Orazio & Duffell (2021). In Sect. 2 we describe the physical assumptions behind the models and the numerical strategy employed to compute the long-term evolution. In Sect. 3 we show the trajectories of individual systems in the eccentricity–mass ratio–semi-major axis parameter space, and identify possible observational signatures for a population of systems. Section 4 is dedicated to the comparison with the models from Zrake et al. (2021) and D’Orazio & Duffell (2021), wherein we highlight key differences and similarities. In Sect. 5 we discuss the implications and applicability of our results to different astrophysical scenarios. Lastly, we summarize our findings in Sect. 6.

2. Methods

We aim to model the long-term orbital evolution of a population of binaries over the lifecycle of disk-mediated accretion phases. Employing hydrodynamic simulations directly for this purpose would be prohibitively expensive. However, we can use the results from fixed-parameter simulations to build a model of the long-term effects of the disk-binary interaction on the orbit. Our methodology involves three steps.

  1. We collect the results of hydrodynamic simulations of circumbinary disks available in the existing literature. From these simulations, we obtain the derivatives describing the rates at which the separation a, eccentricity e, and mass ratio q change as a function of accreted mass (Sects. 2.1 and 2.2, the values used are reported in Appendix C).

  2. We then interpolate these derivatives with a continuous function between the points of the parameter space that were covered by hydrodynamic simulations (Sect. 2.3).

  3. Finally, we use the function obtained in step 2 to numerically solve for the long-term evolution of a system after specifying its initial parameters a0, e0, and q0 (Sect. 2.4).

We release spindler, a python package that performs these calculations, available as a Zenodo release1 and as a github repository2. The code is modular, allowing for future extensions to incorporate findings from additional hydrodynamic simulations, for example, to explore different viscosities or disk aspect ratios. We also release a c-based version of spindler3, for better compatibility with c and Fortran based codes. The interpolation routines of the c version rely on the librinterpolate library4 by Robert Izzard.

2.1. Fiducial suite of hydrodynamics simulations

As a fiducial model, we employ the results of a suite of 2D hydrodynamic simulations carried out and analyzed by Siwek et al. (2023a,b, hereafter S23a and S23b) using the moving-mesh code AREPO (Springel 2010; Pakmor et al. 2016; Muñoz et al. 2013). Here, we summarize the key underlying assumptions.

The simulations represent a finite gas disk accreting onto a central binary on a coplanar prograde orbit of semi-major axis a and eccentricity e. They are performed in 2D, assuming vertical hydrostatic equilibrium with a locally isothermal equation of state. The two objects composing the binary are modeled as sink particles of mass M1 and M2 (where M1 > M2) and total mass M ≡ M1 + M2. These simulations are designed to capture the regime where the objects’ radii are substantially smaller than their separation. Mass is removed from the system from a region within a radius of 0.03a from the sink particles. The sink radius is larger than the objects’ horizon or surface radius, but is determined as small enough to not influence the solutions.

The disk aspect ratio is set to h/r = 0.1, where h is the thickness of the disk and r is the radial coordinate. Viscosity is treated with the α-disk prescription (Shakura & Sunyaev 1973), where Reynolds stresses from unresolved turbulence are characterized by the dimensionless parameter α = 0.1. The gravitational potential of the binary is softened with a softening length of 0.025a, and the mass of the disk is taken to be much smaller than the mass of the binary so that the self-gravity of the disk can be ignored.

Under the stated assumptions, the system of equations becomes scale-free, leaving only two free parameters: the binary mass ratio q ≡ M2/M1 ≤ 1 and the eccentricity of the binary orbit e. S23a explored a rectangular grid of simulations considering ten different mass ratios within the range 0.1 ≤ q ≤ 1, and five values for the orbital eccentricity, spanning from e = 0 to e = 0.8. S23b expanded the grid to nine different eccentricities, spanning the same eccentricity range.

The main quantities we are interested in for this work are the time derivatives of the orbital eccentricity ė, separation ȧ, and mass ratio q ˙ $ \dot{q} $, since these will set the orbital evolution of the binary systems. The binary is not evolved directly in the simulation, as the orbit is kept fixed. Instead, the interaction with the circumbinary disk is tracked by recording the rate of change of angular momentum, energy, and mass at each simulation time-step. From these quantities, the derivatives of the orbital parameters are computed and then averaged over hundreds of orbits, sufficient to remove the effect of the orbital variability, the variability associated with the lump precession (which occurs on a timescale of about five orbits, MacFadyen & Milosavljević 2008; Shi et al. 2012; D’Orazio et al. 2013; Farris et al. 2014; Dittmann & Ryan 2022), and the cavity precession (about 300 orbits, Muñoz & Lithwick 2020; Siwek et al. 2023a).

For convenience, instead of referring to time derivatives, we express the derivatives in terms of the variation of total binary mass:

q ′ ≡ d q d log M , a ′ ≡ d log a d log M , e ′ ≡ d e d log M · $$ \begin{aligned} q{\prime } \equiv {\frac{\mathrm{d}q}{\mathrm{d}\log {M}}}, \quad a{\prime } \equiv {\frac{\mathrm{d}\log {a}}{\mathrm{d}\log {M}}}, \quad e{\prime } \equiv {\frac{\mathrm{d}e}{\mathrm{d}\log {M}}}\cdot \end{aligned} $$(1)

These expressions have the advantage of being invariant under the rescaling of the simulation’s physical units of binary mass and gas surface density.

Our approach to compute the long-term orbital evolution is valid when in a quasi-equilibrium state, that is the limit where the parameters of the binary vary on a timescale much longer than the viscous relaxation time tvisc of the inner region of the disk, which dominates the forces on the binary. This condition is enforced in the simulations, but it depends on the astrophysical system of interest whether it is satisfied.

2.2. Additional simulations

To better understand the impact of model uncertainties, we compare our main results, based on hydrodynamic simulations by S23b, with alternative simulations by Zrake et al. (2021, hereafter Z21) and D’Orazio & Duffell (2021, hereafter DD21). These works modeled circumbinary accretion and the effect on the orbit under nearly identical assumptions. The main difference from S23b is that these two studies focus exclusively on q = 1 binaries, but with a far denser sampling in eccentricity. Below we detail the other modeling variations, which are additionally summarized in Table 1.

Table 1.

Summary of the main parameters and results of the considered suites of hydrodynamic simulations.

Beyond the differences in parameter exploration, the primary distinction between the three works is the structure of the underlying simulation grid. S23b uses the Voronoi moving mesh of Arepo, Z21 employs the Cartesian grid-based code with static mesh-refinement Mara3 (Zrake & MacFadyen 2012; Tiede et al. 2020), and DD21 rely on the cylindrical moving mesh of DISCO (Duffell 2016; Duffell et al. 2020). Each technique has its respective strengths. Moreover, compiling the results bolsters confidence in features on which the three agree. Conversely, our approach diminishes the likelihood of capturing dynamics that may be emergent from the underlying mesh structure.

The treatment of viscosity in the three studies is nearly identical. S23b and Z21 both employ an α-disk prescription (Shakura & Sunyaev 1973), setting α = 0.1 (see the Appendix of S23a on the implementation for a binary potential). DD21 assumes a constant kinematic viscosity ν = 10−3a2Ω, where Ω = 2π/Porb is the orbital frequency and Porb is the orbital period. The two prescriptions yield a comparable value of viscosity at r = a, for h/r = 0.1.

All three studies treat accretion onto the individual components of the binary with a sink particle at the location of the binary objects. The treatment of the sink particles can affect the morphology of the “mini-disks” formed around them, as well as the rate at which mass and angular momentum are imparted onto the binary. All three studies employ the commonly used sink prescription for which the removed angular momentum is the one associated with the removed mass at the local fluid velocity (as opposed to “torque-free” or “spinless” sinks, Dempsey et al. 2020). As such, all of these studies will possess some induced spin-component in the binary angular momentum (and this spin component would be real in the limit of a resolved stellar surface or black hole Innermost Stable Circular Orbit, ISCO). This could be relevant for some subsets of accreting binaries, and while the effects of this spin component have been documented for circular binaries (e.g., Muñoz et al. 2019; Dittmann & Ryan 2021), they have not been examined in detail for eccentric orbits. The sink radius is always taken as a small fraction of the binary semi-major axis (see Table 1).

The studies also differ in the strategy for exploring the parameter space and in the duration of the simulations. Z21 and S23b run a discrete number of individual simulations. Each simulation concerns a particular binary eccentricity and mass ratio, that are kept fixed for the whole run. Each simulation lasts 2000 orbits and 10 000 orbits, respectively. At the end of the run, the forces acting on the binary and other quantities of interest are computed and averaged over the last hundreds of orbits. Z21 exclusively explore the q = 1 regime, but do so with 32 simulations and a fine eccentricity spacing. In comparison, S23b conduct nine simulations at q = 1.

DD21 took a distinct approach, running a single extended simulation lasting 20 000 orbits. In this simulation, eccentricity is gradually increased from e = 0 to e = 0.9. In this work we only consider the part with e ≤ 0.8 in order to be consistent with the samples taken from the other two studies, and because for e > 0.8 the size of the mini-disks at the periapsis passage becomes smaller than the gravitational smoothing length, as noted by the authors. The resulting a′ and e′ are then smoothed over an ≈450-orbits window to average over the short-term variability and produce a result equivalent to having run ≈40 simulations with eccentricity spaced by Δe = 0.02. As a consistency check, DD21 have run additional simulations at a fixed orbital configuration, and compared them with the values obtained through the eccentricity sweep, finding a negligible difference in most cases. The advantage of the approach by DD21 is that, since the eccentricity is increased slowly, the disk is kept close to a relaxed state at all times, preventing the need to spend computational resources in relaxing a new disk for each eccentricity.

2.3. Interpolation of the derivatives

To compute the long-term orbital evolution of the binary systems, we linearly interpolate the derivatives a′, e′ and q′ (as defined in Eq. (1)) obtained from the hydrodynamic models. To this aim, we define a set of interpolating functions fa(e, q), fe(e, q), and fq(e, q). These functions do not depend on the separation a because of the scale invariance of the problem. For our fiducial model, based on S23b, we use bilinear interpolation5 (Weiser & Zarantonello 1988) on the rectangular grid of eccentricities and mass ratios to obtain the functions fa(e, q) and fe(e, q), which interpolate the values provided in Tables C.1 and C.2. These are the same values reported in Figs. 2 and 3 of S23b, with an additional suite of simulations at e = 0.05, previously unpublished.

We obtain fq(e, q) first by interpolating the relative accretion rates λ ≡ Ṁ2/Ṁ1 published in Fig. 4 of S23a (also provided in our Appendix in Table C.3 for convenience), then by computing

f q ( e , q ) ≡ ( 1 + q ) ( λ ( e , q ) − q ) 1 + λ ( e , q ) · $$ \begin{aligned} f_q(e,q) \equiv \frac{(1+q)(\lambda (e,q) - q)}{1+\lambda (e,q)}\cdot \end{aligned} $$(2)

Our functions fa, fe and fq are defined for 0 ≤ e ≤ 0.8 and 0.1 ≤ q ≤ 1. In case the integration routine requires to evaluate them outside these ranges, we return the value of the nearest known data point.

We tested a cubic spline and the PCHIP method (Fritsch & Butland 1984) as alternative interpolation techniques. We observed that using a cubic spline introduced spurious oscillations in the interpolating functions. Instead, the PCHIP method, which is guaranteed to preserve the monotonicity of the data, gave satisfactory results that closely resembled linear interpolation. Since PCHIP and the linear method gave equivalent results, we opted to use the simpler linear interpolation.

For Z21 and DD21 we only have to consider one dimension: eccentricity. In the case of Z21, the derivative of the semi-major axis is not provided in their work. We publish it here, together with the eccentricity derivative, and provide them in Table C.4, in Appendix C. In the case of DD21, their eccentricity sweep yielded continuous curves for the derivatives a′ and e′ instead of discrete points. Their results show oscillations (see Fig. 1 in DD21), even after their smoothing procedure. As we expect that these oscillations do not represent the steady state, we remove them by discrete resampling with carefully selected points to preserve the major features as shown in Fig. B.1. The results are reported in Table C.4. We then obtain our interpolating functions for Z21 and DD21 by linear interpolation of these tables.

2.4. Integration of the orbital evolution

The evolution of the binary system is then computed by solving the differential equations

a ′ ≡ d log a d log M = f a ( e , q ) , $$ \begin{aligned}&a{\prime } \equiv {\frac{\mathrm{d}\log a}{\mathrm{d}\log M}} = f_a(e,q),\end{aligned} $$(3)

e ′ ≡ d e d log M = f e ( e , q ) , $$ \begin{aligned}&e{\prime } \equiv {\frac{\mathrm{d}e}{\mathrm{d}\log M}} = f_e(e,q),\end{aligned} $$(4)

q ′ ≡ d q d log M = f q ( e , q ) , $$ \begin{aligned}&q{\prime } \equiv {\frac{\mathrm{d}q}{\mathrm{d}\log M}} = f_q(e,q), \end{aligned} $$(5)

for the variables a, e, and q, with appropriate initial conditions e0 and q0. Given the scale invariance of the problem, we always set a0 = 1. The integration was performed with the LSODA method (Petzold 1983) as implemented in SciPy. The relative tolerance (rtol) was set at 10−2, and the absolute tolerance (atol) at 10−5.

3. Results

3.1. Individual binary evolution

The evolution of a binary interacting with and accreting from a circumbinary disk is captured by the derivatives a′, q′ and e′ (Eqs. (3)–(5)), which describe the change of each of these parameters as a function of the accreted mass.

In Fig. 1, we show our interpolation for the evolution of the semi-major axis as a function of eccentricity and mass ratio. The red regions denote values of q and e where the disk increases the binary separation, while the blue regions mark where it causes the binary separation to shrink. The white lines with arrows show the trajectories traced by example systems as they evolve in q-e space. The black band in the plot center denotes the equilibrium eccentricity as a function of q (see S23b). All paths first move toward the equilibrium eccentricity, after which they tend to grow their mass ratio toward unity. In the limit of infinite mass accretion, all streamlines eventually converge on a single attractor, with a mass ratio qatt = 1 and eccentricity eatt ≈ 0.5. At the attractor, binaries tend to shrink their orbit at a rate a att ′ =−0.7 $ a^\prime_{\rm att} = -0.7 $, or equivalently, like a ∝ M−0.7.

thumbnail Fig. 1.

Orbital evolution of binary systems under interaction with a circumbinary disk. The white lines follow the trajectories in the eccentricity–mass ratio plane. The thick black line traces the equilibrium eccentricity, that is the region where, for a fixed mass ratio q, the eccentricity derivative is e′ = 0. The background color shows our interpolated function approximating the derivative of the separation a′≡dlog a/dlog M. This quantity shows how fast the separation changes per unit of accreted mass. Regions where the orbit shrinks are shown in blue and regions where the orbit widens in red. The small black dots in the background correspond to the parameters used in the hydrodynamic simulations by S23b, on which the interpolation is based.

In a realistic scenario, systems may not reach the attractor because the mass reservoir feeding the disk is depleted or because other mechanisms begin to dominate the evolution of the orbit, such as stellar tides or gravitational-wave emission. In Fig. 2 we quantify how much mass a system needs to accrete in order to approach both the equilibrium eccentricity and the final attractor, in the bottom and top panels, respectively. Systems typically reach the equilibrium band after accreting less than 10% of their initial mass, while they may require up to three times the initial mass to reach the attractor point at q = 1.

thumbnail Fig. 2.

Effect of accretion on the orbit of the binary. The top panel shows how much mass the systems need to accrete to reach the vicinity of the attractor point (white circle). The bottom panel shows the amount of accreted mass needed to reach the equilibrium band only (white region). The black line traces the equilibrium eccentricity. See also Fig. 1.

Binaries tend to adjust their eccentricity toward the equilibrium value before meaningfully altering their mass, as shown by the nearly vertical arrows in Fig. 1. This is because the magnitude of e′ is several times larger than q′. This is shown in more detail for the same streamlines in Fig. 3 where the coloring shows the evolution rates for q and e in the top and bottom panels respectively. It is important to note the difference in the dynamic range (max|e′| ∼ 9 while max|q′| ∼ 1.5), which can be read from the color bar. Figure 3 shows that the systems whose eccentricity evolves faster are those with nearly equal masses and low eccentricities (0.8 < q < 1, 0.2 < e < 0.4), and those with large eccentricities and mass ratios (q < 0.2, e > 0.6). The least sensitive systems are nearly circular (e = 0 − 0.1) and those near the equilibrium eccentricity band. There also exists an asymmetry for the behavior of systems with near equal masses (q > 0.8) in that the pumping of e at initially modest eccentricities (0.1 < e < 0.4) is much stronger than the damping in systems with initially high e ≳ 0.5 (see also the discussions in Dunhill et al. 2015 and S23a).

thumbnail Fig. 3.

Velocity at which systems travel across the eccentricity–mass ratio plane. The streamlines are the same as in Fig. 1, but here the color of the streamlines encodes how fast the mass ratio (top) and eccentricity (bottom) change per unit of accreted mass, quantified by the derivatives |q′| and |e′|. The black line traces the equilibrium eccentricity.

Generally, only when binaries have functionally reached the equilibrium eccentricity (or are almost perfectly circular) does the evolution of q start to dominate, |q′| > |e′|. As a result, the evolution of a typical system is composed of two phases. In the first phase, eccentricity evolves toward the equilibrium band, keeping the mass ratio almost constant. In the second phase, the mass ratio moves toward unity, while the eccentricity tracks the q-dependent equilibrium value.

Figure 4 shows the evolution of the eccentricity and semi-major axis of a series of accreting binary systems. Each system is initialized with a different eccentricity along the vertical gray bar on the right side of each panel, and the initial mass ratio is noted in each panel’s upper left corner. As mass is accreted on to the binary, the eccentricity again changes most rapidly such that each system initially reaches the equilibrium value after accreting ≈10% of its original mass, and then subsequently only varies with the q-dependence of the equilibrium eccentricity as the binary mass ratio grows.

thumbnail Fig. 4.

Evolution of eccentricity (e) and semi-major axis (a) as a function of accreted mass. For each panel, the systems start from a different initial mass ratio (q0) that spans from 0.1 (top left) to 1.0 (bottom right). Within each panel, the initial eccentricity (e0) spans from 0.01 to 0.8. The trajectories are colored based on the mass change ΔM/M0, where M0 is the initial mass of the binary and ΔM is the amount of accreted mass. All the systems start along the vertical line a = a0, highlighted with a thick gray line.

Despite the presence of several regions of the e-q space with a positive a′, in the region of equilibrium eccentricity a′ is almost always negative, and the orbit will tend to shrink. However, similar to the mass ratio case, in order to achieve a significant reduction of the orbital separation |Δa|≈a0, the system needs to accrete from the disk an amount of mass comparable to the initial mass of the system ΔM ≈ M0. This implies that orbital decay from a steady, disk-mediated accretion phase may only be viable in the limit that the gas reservoir is at least of comparable mass to the binary itself.

3.2. Population evolution

In addition to examining the evolutionary trajectories of individual accreting systems, we also explore the properties of populations and describe the observational signatures that they produce on the distributions of orbital parameters. We randomly sample a population of binary systems from a flat distribution of eccentricity between 0 ≤ e0 ≤ 0.8 and a flat distribution of mass ratio between 0.1 ≤ q0 ≤ 1. Because the problem is scale-free, all systems are initialized with a semi-major axis of a0 = 1. Then, assuming each system accretes at the same rate, we model our population forward in time according to the solutions of Eqs. (3)–(5). We choose these simple initial conditions to be agnostic to the type of binary systems (stars, black holes, etc.) that we are modeling. We sample a population of 106 systems, which we deem large enough to densely cover the two-dimensional parameter space of eccentricity and mass ratio.

Figure 5 shows the population distributions of eccentricity, semi-major axis, and mass ratio at different accreted masses ΔM. Up to the accretion of ≈10% of the initial mass (ΔM/M0 ≤ 0.1), the most noticeable variation happens in the eccentricity distribution in the top panel, which concentrates around e ≈ 0.5 After accreting more than 50% of the initial mass, however, the mass ratio and semi-major axis distributions begin to change noticeably. The mass ratio distribution develops a peak at q = 1 as the initially largest-q systems grow to equal mass, and mass ratios below q = 0.5 are depopulated. The distribution of a/a0 broadens to ≈0.6 dex and features a double peak. The rightmost peak is composed of those systems that start with q > 0.8, which experience a slow evolution of the semi-major axis due to a small value of a′(e  =  eeq)≈ − 0.7. On the contrary, the systems that start with q < 0.8 can have a′(e  =  eeq) as large as −4 (see also Fig. B.2 in Appendix B) so they migrate to a tighter orbit faster, and eventually accumulate on the leftmost peak upon reaching equal mass ratio.

thumbnail Fig. 5.

Evolution of the distribution of the orbital parameters for a population of binaries. Distribution of eccentricity (top), semi-major axis (middle), and mass ratio (bottom). Each color corresponds to a different amount of accreted mass, increasing from dark purple to light yellow.

Figure 6 depicts the evolution of the population in the e-q plane. In the top-left panel, every system has accreted 1% of the initial mass, and they remain close to their initially flat distributions because 1% is not sufficient to alter the distribution of eccentricity or mass ratio. After accreting 5% of the initial mass – in the top-right panel – a pattern starts to emerge and becomes evident in the middle-left panel as a clear correlation between eccentricity and mass ratio, caused by the accumulation of systems close to their equilibrium eccentricity. When the systems have accreted 20% of the initial mass (middle-right panel) the e-q correlation becomes extremely tight. In the remaining two panels, after accreting 50% or 100% of the initial mass, all the systems settle at the same eccentricity e ≈ 0.5 and move toward a mass ratio of one. This e-q correlation could be observed in populations of binaries that have undergone steady, disk-mediated accretion in the past, as long as additional dynamics have not since meaningfully modified the orbit. The presence of the correlation could also help constrain the typical amount of mass that was accreted from the disk.

thumbnail Fig. 6.

Population synthesis of binaries interacting with circumbinary disks in the e-q plane. Each panel corresponds to a different amount of accreted mass, increasing from top-left to bottom-right. The color scale indicates the number of systems falling into each bin. The gray band marks the position of the equilibrium eccentricity.

Because of the long real-time evolutionary timescales of the relevant stellar or black hole binaries, identifying (or ruling out) features in parameter distributions of populations is likely the most promising approach for determining and constraining the importance of accretion phases in binary evolution.

4. Comparison with other models

Here we compare the results of our fiducial model to those based on the additional simulations discussed in Sect. 2.2. The predictions of the three models broadly agree, all of them developing an equilibrium eccentricity around e ≈ 0.4 − 0.5. An intuitive explanation for the existence of an equilibrium eccentricity can be understood in terms of the interaction between the inner edge of the disk and the less massive component of the binary, as described by Artymowicz et al. (1991) and Roedig et al. (2011). Namely, when the binary is at apoapsis, the less massive object approaches the inner disk edge. If its angular velocity surpasses that of the fluid in the disk, it is decelerated, leading to an increase in orbital eccentricity. Conversely, when the object’s angular velocity lags behind that of the disk fluid, the object is propelled forward, subsequently reducing the orbital eccentricity. Equilibrium is approximately achieved when the angular velocities of the object at apoapsis and the disk fluid are the same. We note, however, that this intuitive picture is approximate and does not fully account for the disk’s morphological complexities in the full nonlinear solutions. For instance, the equilibrium eccentricity in DD21 does not conform to this picture, but is associated with the disk transitioning from a symmetric to an asymmetric state.

We report here (and in Table 1) the equilibrium eccentricity for the three models at q = 1, obtained via linear interpolation. We find e eq = 0 . 48 − 0.08 + 0.02 $ e_{\mathrm{eq}} = 0.48^{+0.02}_{-0.08} $ for S23b, 0 . 443 − 0.018 + 0.007 $ 0.443^{+0.007}_{-0.018} $ for Z21, and 0.39 ± 0.01 for DD21. The uncertainty is taken to be the distance from the closest known point, except in the case of DD21, where it is half of the reported eccentricity resolution Δe = 0.02. These results are consistent with the ones reported in the respective papers. See also Fig. B.3 for a comparison of the derivatives obtained by the three suites of simulations.

4.1. Behavior at small eccentricity

The primary evolutionary difference between the three suites of simulations is the sign of e′ for e ≲ 0.1. Z21 and DD21 both determine that such near-circular orbits are stable and remain near-circular with e′< 0, while S23b find that such near-circular orbits are unstable with e′> 0 always. This small discrepancy implies a dramatically different evolution for systems starting from nearly circular orbits. According to S23b, those systems will increase their eccentricity up to the equilibrium value e ≈ 0.5 and then shrink the orbit, while for Z21 and DD21, systems with small eccentricity will tend to circularize while they grow their orbital separation. A positive e′ for e ≤ 0.2 was also found by SPH studies (e.g., Artymowicz et al. 1991; Cuadra et al. 2009; Roedig et al. 2011; Fleming & Quinn 2017).

The reason for this discrepancy is currently unclear. S23b attributes this discrepancy to different radial dependencies in the kinematic viscosity adopted, for instance by DD21, but this cannot explain how S23b and Z21 – who adopt the same viscosity prescription – obtain a different result. Aside from the grid structure, the other major difference between these studies is that the mass removal rate in the sink-region (see Sect. 2.2) is much more rapid in S23b than in the two additional studies. However, some of us have examined whether this could account for the discrepancy and found that commensurately increasing the sink-rate in Eulerian grid-based simulations like DD21 and Z21 actually exacerbates the disparity. See Sect. 5.3 for a discussion of the consequences for post–common-envelope systems.

4.2. Orbital decay

Figure 7 compares the long-term evolution of individual systems according to the three models, in the regime of q = 1. The top row illustrates the change in eccentricity. Each system starts from a different initial eccentricity e0, spanning from 0.01 to 0.8, encoded by color. The three models agree that after accreting ≈10% of the initial mass (ΔM/M0 ≈ 0.1), all systems converge toward an equilibrium eccentricity. The only exception is the system initialized with e0 = 0.01 (dark blue line), which tends to circularize for Z21 and DD21, since e′< 0 in those models.

thumbnail Fig. 7.

Model comparison of the orbital evolution of equal-mass binaries. We show eccentricity (top) and semi-major axis (bottom) as a function of accreted mass. Each line represents a different system, colored according to the initial eccentricity, in the range 0.01 ≤ e0 ≤ 0.8. The black dashed lines in the bottom panels mark the analytical solutions obtained in Sect. 4.2 considering a system evolving at the equilibrium eccentricity.

The bottom row shows the evolution of the semi-major axis. The left panel confirms that, according to S23b, all the systems are led to orbital shrinking independently of the initial eccentricity. Conversely, the central and right panels show that for Z21 and DD21, the orbit shrinks only if the initial eccentricity is larger than about 0.1. The amount of accreted mass necessary to change the orbital separation is also model-dependent. This is due to the different values of a′ at the equilibrium eccentricity in the three models (see Table 1). The black dashed lines correspond to the case of a system evolving at the equilibrium eccentricity, whose evolution is given by the solution to the equation

d log a d log M = γ ( constant ) ⇒ a ∝ M γ , $$ \begin{aligned} {\frac{\mathrm{d}\log a}{\mathrm{d}\log M}} = \gamma \quad \mathrm{(constant)} \implies a \propto M^\gamma , \end{aligned} $$

where γ is the value of the derivative of the semi-major axis at the equilibrium eccentricity. The linear interpolation performed on the results of S23b and Z21 gives a′(e  =  eeq) = − 0.70 and a′(e  =  eeq) = − 2.3, respectively. Conversely, DD21 report a′(e  =  eeq)≈0, which the authors explain with the presence of an abrupt transition in the disk geometry in correspondence with the equilibrium eccentricity. However, DD21 predict oscillations around the equilibrium configurations to yield an effective a′(e  =  eeq) = − 0.49 (see also Fig. B.2). We note that the present work does not take into account such oscillations when computing the orbital evolution, hence the mismatch between analytical and numerical prediction in the bottom right panel of Fig. 7.

The different values of a′(e  =  eeq) obtained by the three models imply a different amount of orbital shrinking, for a given accreted mass. Namely, after doubling the mass of the system, we predict the semi-major axis to have shrunk by a factor 5, 1.6, and 1.4, for Z21, S23b and DD21, respectively.

4.3. Population comparison

Figure 8 compares the distribution of the orbital parameters of a population of binary systems interacting with circumbinary disks, according to the three models. Here, similarly to Sect. 3.2, we have initialized 106 systems with a flat distribution of eccentricity from 0 to 0.8 and a constant mass ratio q = 1. The top row presents the distribution of eccentricity. For S23b, all the systems transition toward the equilibrium configuration, producing a visible peak in the eccentricity distribution already after the accretion of 2% of the initial mass. The middle panel, relative to Z21, highlights that the eccentricity evolves toward a bimodal distribution peaked at e = 0 and e ≈ 0.5, with a gap at e ≈ 0.1 (see also Fig. 1 of Z21). The right panel shows the result for DD21, which predicts an additional third peak, due to an accumulation of systems around e ≈ 0.2. However, the third peak disappears when ΔM/M0 ≥ 0.1.

thumbnail Fig. 8.

Model comparison of the evolution of the distributions of the orbital parameters for a population of equal-mass binaries. Distribution of eccentricity (top), semi-major axis (bottom), colored based on the amount of accreted mass. Each column shows one of the three models.

The bottom panels show the distribution of the semi-major axis divided by the initial value a0. All the systems start with a/a0 = 1, but after having accreted some mass from the circumbinary disk, the semi-major axis distribution spreads. For S23b, all the systems are led to orbital shrinking, and this is reflected in the distribution shifting toward shorter semi-major axes. Z21 presents a similar trend, except for the ≈10% of the systems that started with e0 < 0.1, so that their orbit is widening. The model based on DD21 also predicts the systems with small eccentricity to widen, while the others maintain approximately a constant semi-major axis. As discussed for Fig. 7, when considering small oscillations around the equilibrium eccentricity, the DD21 model would instead predict orbital shrinking.

5. Astrophysical implications

We have investigated the long-term orbital evolution of binary systems subjected to the accretion of a circumbinary gas disk by employing the results of three suites of 2D hydrodynamic simulations. The scale-free nature of these simulations encompass a diverse range of astrophysical scenarios, including binary star formation (Sect. 5.1), interacting binaries (Sect. 5.2), post–common-envelope systems (Sect. 5.3), post-AGB binaries (Sect. 5.4) and SMBH (Sect. 5.5). We discuss here the astrophysical and observational implications, and describe how these results can help guide the direction of future efforts and advancements in the field.

5.1. Binary star formation

One proposed mechanism for the formation of close binary stars (see Offner et al. 2023, for a review) revolves around the concept of disk fragmentation (Adams et al. 1989; Shu et al. 1990; Bonnell & Bate 1994; Mignon-Risse et al. 2023). In this scenario, gravitationally unstable clumps in the proto-stellar disk undergo inward migration as they grow (Ward 1997; Kley & Nelson 2012), ultimately resulting in the formation of a central gas cavity around the two young stellar objects. After the formation of the cavity, if the orbital evolution is primarily governed by the interaction with a stable thin accretion disk, similar to the ones considered in this work, we expect to observe clear signs of this interaction in the eccentricity and mass ratio distributions among young main-sequence binaries. Particularly, we expect this signature to be evident in systems with sufficiently wide orbits, which have remained unaffected by tides. According to our fiducial model, we expect a correlation between eccentricity and mass ratio, as illustrated in Fig. 6. In contrast, the models based on Z21 and DD21 predict a bimodal eccentricity distribution with peaks at e = 0 and e ≈ 0.4, corresponding to the two equilibria.

Observationally, such Young Stellar Objects (YSOs) come with a variety of different disks – with varying adherence to the assumptions laid out in Sect. 2 – the most pertinent aspects of which we detail below and summarize in Table 9. YSOs can be subdivided into classes (Dunham et al. 2014) depending on the shape of their spectra, which also roughly correspond to progressive evolutionary states.

thumbnail Fig. 9.

Summary of the properties of circumbinary disks in astrophysical scenarios. Binary mass refers to the initial mass M0 of the central binary. The disk is fed if there is a mass inflow that replenishes the disk. The disk is self-gravitating if the self-gravity of the disk has a nonnegligible effect on the dynamics and structure of the disk. Disk mass Md is the instantaneous mass of the disk. Accretion rate is the mass accreted by the binary per unit time. Disk lifetime is the duration of a disk-mediated phase in the binary evolution, terminating when the disk is dispersed or the binary merges. Accreted mass is the amount of mass that can be accreted from the disk on to the binary within the disk lifetime. When the disk is not fed from external mass inflow, the accreted mass can be estimated as Md/M0, which assumes that all the mass in the disk is eventually accreted, and neglects outflows or other mechanisms of mass loss. h/r is the disk aspect ratio in the vicinity of the binary. α is the viscosity parameter of the disk. [1]Jørgensen et al. (2009), [2]Tobin et al. (2020), [3]Jacobsen et al. (2018), [4]Kido et al. (2023), [5]Fiorellino et al. (2023), [6]Takakuwa et al. (2014), [7]Manara et al. (2020), [8]Guilloteau et al. (1999), [9]Scaife (2013), [10]Dutrey et al. (2016), [11]Long et al. (2021), [12]Testi et al. (2022), [13]Grant et al. (2023), [14]Tuna & Metzger (2023), [15]Vos et al. (2015), [16]Izzard & Jermyn (2023), [17]Bujarrabal et al. (2018), [18]Kluska et al. (2022), [19]Bujarrabal et al. (2013), [20]Franchini et al. (2021), [21]Novak et al. (2011), [22]Schawinski et al. (2015), [23]Siwek et al. (2020), [24] estimated from Eq. (A.1) assuming r = 100 AU and α = 0.1, [25] estimated from Eq. (A.2) assuming r = 0.01 pc, α = 0.1 and η = 0.06, [26]Rosotti (2023), [27]Starling et al. (2004).

Class 0 YSOs are still too deeply embedded in gas to be detectable in the near infrared. Measuring their properties is challenging, but it is often possible to determine the presence of an accretion disk with a mass around 1–10% of the envelope mass (Jørgensen et al. 2009). In this regime, the disk is expected to be thick, self-gravitating, and possibly unstable to fragmentation (Kratter & Lodato 2016). For example, the inward migration and subsequent accretion of fragments by the central star has been proposed as the mechanism driving the outbursts of FU Orionis stars (Hartmann & Kenyon 1996; Elbakyan et al. 2021). This merger-dominated mass assembly has also been invoked as a likely pathway to the formation of binary stars (Tokovinin & Moe 2020; Riaz et al. 2021). Therefore, in the early phases of star formation, when accretion is bursty and dominated by mergers, the steady accretion models employed in this work may not be applicable.

In Class I YSOs (the next evolutionary stage), the envelope has largely dissipated and 70–98% of the mass in the star-disk-envelope system is contained in the star itself (Jørgensen et al. 2009). Therefore, assuming no mass inflows or outflows, we estimate a limit on the amount of mass that can still be accreted during star formation to ΔM/M0 ≈ 2 − 40%. Fiorellino et al. (2023) find that 20–30% of the YSOs observed in this stage have a disk mass too small to trigger fragmentation, resulting in stable disks like the ones considered in this work. One of the few examples of a binary system observed in this configuration is L1551 (Takakuwa et al. 2014).

Class II YSOs and T Tauri stars correspond to the next evolutionary stage. At this stage, the envelope has completely dissipated and it is usually assumed that there is no significant inflow to feed the disk (although recent findings are challenging this assumption, see Gupta et al. 2023; Kuffmeier et al. 2023). As a consequence, the mass accreted in this stage is limited by the mass of the disk, which is typically 10−4 − 10−1 M⊙ (e.g., Manara et al. 2020; Testi et al. 2022). Assuming that all the matter in the disk is eventually accreted by the star gives ΔM/M0 ≈ 10−3 − 10−1. In practice, the accreted mass will be smaller, due to photoevaporation (Somigliana et al. 2020; Winter & Haworth 2022) or magnetized winds (Lesur 2021; Turpin & Nelson 2024). Arguably, the most studied Class II binary system interacting with a circumbinary disk is GG Tau (Dutrey et al. 1994; Guilloteau et al. 1999; Scaife 2013; Andrews et al. 2014; Aly et al. 2018), but see also GM Au (Dutrey et al. 1998), TWA 3A (Tofflemire et al. 2017) the Herbig stars HD 142527 (Lacour et al. 2016; Price et al. 2018) and V892 Tau (Long et al. 2021).

Having revised the disk properties of YSOs, we conclude the following: In the earliest phases of star formation, circumbinary disks are thick, self-gravitating or fragmenting, and result in high accretion rates (> 10−6 M⊙ yr−1) onto their YSOs. They later transition to a stable, thin configuration with a much lower accretion rate (< 10−7 M⊙ yr−1). The amount of mass assembled in each phase, and the relative importance of steady versus bursty accretion is not yet completely understood. However, as long as the final few percent of the mass is accreted from a stable thin circumbinary disk, during the Class I or Class II phases, it should be possible to discern the signatures of an interaction with the disk predicted in this work. To our knowledge, such features in the main-sequence star populations have not been reported in the literature (e.g., Tokovinin 2020; Hwang et al. 2022; Andrew et al. 2022; Gaia Collaboration 2023).

However, there exist some plausible explanations for the lack of observational evidence for our predicted signatures. In the first place, the mass accreted from a thin, stable disk during star formation could account for less than ∼1% of the mass of the binary. In this case, the effect of the interaction with the disk would be negligible (see Figs. 6 and 7). In the second place, some physical processes currently not included in our model (such as outflows, radiation transport, or dynamical interactions with other stars) may play an important role in the orbital evolution of binary YSOs. In the third place, the evolution might diverge from our predictions if the assumption of a relaxed disk becomes invalid (see Sect. 2.1). Our approach assumes that the orbital parameters evolve on a longer timescale than the viscous time tvisc(r) = 2/3 r2/ν(r) (e.g., Haiman et al. 2009). For the inner region r = 2a, which dominates the forces, tvisc ≈ 300Porb. We estimate the characteristic timescale for the change of orbital parameters as

t X ≡ 1 X ′ M M ˙ = 1000 P orb 1 X ′ ( a 100 AU ) − 3 / 2 ( M M ⊙ ) 3 / 2 ( M ˙ 10 − 6 M ⊙ yr − 1 ) − 1 , $$ \begin{aligned} t_X &\equiv \frac{1}{X{\prime }} \frac{M}{\dot{M}} = 1000 P_{\mathrm{orb}} \frac{1}{X{\prime }} \left(\frac{a}{{100}\,\mathrm{AU}}\right)^{-3/2}\nonumber \\& \quad \left(\frac{M}{{M}_{\odot }}\right)^{3/2} \left(\frac{\dot{M}}{10^{-6}\,{M}_{\odot }\,\mathrm{yr}^{-1}}\right)^{-1}, \end{aligned} $$(6)

where X stands for a, e or q and X′ is defined as in Eq. (1) and can be as high as ∼8 (see Fig. 3). Hence, in systems with a wide orbit and a high accretion rate, tX can be comparable or even shorter than tvisc. This violates the assumption that the disk is relaxed, especially for wide binaries, where tidal effects are negligible. Furthermore, the hydrodynamic simulations used herein fix the binary orbit and do not account for the back-reaction from the surrounding gas (see e.g., Franchini et al. 2023; Tiede et al. 2024, for discussion) In the fourth place, the assumption that the size of the binary components is much smaller than the orbit may become invalid for short-period systems and when the effective accretion radius of the stars becomes large (e.g., because of the presence of a magnetosphere). A fifth caveat concerns the treatment of viscosity. Magnetorotational instability, which is generally believed to source the Reynolds stresses that redistribute angular momentum in many disks, is significantly suppressed in cool proto-stellar disks (King et al. 2007). Measurements of viscosity for Class II protostars indicate 3 × 10−4 ≤ α ≤ 3 × 10−3 (Rosotti 2023), much lower than the value α = 0.1 employed in this work. A final possibility is that only a fraction of the binary systems undergo stable, thin, disk accretion. In this case, the observational signatures would not immediately stand out in the statistics of a general population. A similar idea was brought forward by El-Badry et al. (2019), who proposed that if a subset of the binaries accreted a significant fraction of their mass from a circumbinary disk, they would appear as a population of equal mass “twin” binaries and explain the excess of q ≈ 1 systems in the observed distributions (see also the discussion in Duffell et al. 2020).

Before concluding the discussion about the astrophysical implications of our models to binary star formation, we also mention that binary-disk interactions have been hypothesized as a possible mechanism to explain the observed increase of radial velocity dispersion in young clusters as they age (Sana et al. 2017; Ramírez-Tannus et al. 2017, 2020, 2021). However, to explain the considerable decrease in the minimum period, from 3500 d in the young cluster M17 to 1.4 d in the older clusters, a binary with an initial total mass of M = 20 M⊙ should decrease their separation by more than two orders of magnitude. Based on our current findings, achieving such a significant reduction in separation from a stable accretion scenario would require the system to increase its mass nearly tenfold, while the mass available in the disk is typically less than 10−3 M⊙ (e.g., Backs et al. 2023).

5.2. Mass transferring stars

Circumbinary disks can also be formed during nonconservative mass transfer episodes. These circumbinary disks have been studied in the context of cataclysmic variables (Spruit & Taam 2001; Taam & Spruit 2001; Dubus et al. 2002) and they have been observed around X-ray binaries (Blundell et al. 2008). We refrain from giving quantitative statements regarding the orbital evolution in these cases, since the gas does not accrete on the central binary, but emanates from it. Given this crucial difference, our models are not directly applicable. However, this scenario is worth studying with detailed hydrodynamic models, as it has important consequences on the evolution of the binary and, for example, on the merger rate of gravitational-wave sources (e.g., van Son et al. 2022).

A further occurrence of circumbinary disks is found during stable mass transfer in a triple stellar system, where in some cases the mass stream from the outer star can form a circumbinary disk around the inner binary (de Vries et al. 2014; Leigh et al. 2020; Di Stefano 2020; Dorozsmai et al. 2024). There is currently no direct observational counterpart of this process, but it is expected to be a natural stage of triple star evolution (Tauris & van den Heuvel 2014; Hamers et al. 2022). The effect of the disk on the orbit during triple mass transfer has often been neglected in past population synthesis studies. We believe that spindler will help to improve on this aspect, by providing a simple interface to access the angular momentum and orbital energy exchanged with the disk.

5.3. Post–common-envelope systems

A scenario wherein the interaction with a circumbinary disk may play a crucial role is during the post–common-envelope phase. Following the initial dynamical phase of plunge-in, 3D simulations show that the inspiral process may stall even if the envelope has not been completely ejected. The stalling is caused by the formation of a low gas-density region at the envelope’s center, where friction is reduced (Passy et al. 2012; Ricker & Taam 2012; Ricker et al. 2019; Moreno et al. 2022; Lau et al. 2022a,b). Subsequently, the bound part of the remaining envelope is expected to cool and form a centrifugally supported disk, although the process takes too long to be captured in hydrodynamic simulations of the post–common-envelope phase (Reichardt et al. 2019; Gagnier & Pejcha 2023, 2024).

The interaction of the binary with the disk has been proposed as a mechanism to further reduce the orbital separation (Kashi & Soker 2011; Tuna & Metzger 2023; Wei et al. 2024), which could potentially lead to a late merger. Consequently, the final state of the system could be meaningfully influenced by the interaction with a disk, carrying crucial implications for double compact object mergers and gravitational-wave detection rates (e.g., Belczynski et al. 2007; Stevenson et al. 2017).

Although the eccentricity of the binary before the common envelope can be arbitrarily high (e.g., Vigna-Gómez et al. 2020), it is expected that when the plunge-in stalls the orbit will have become nearly circular, with only a small residual eccentricity e < 0.2 (Trani et al. 2022). Indeed, observed post–common-envelope candidates (identified as short period binaries with at least one hydrogen-depleted component) appear to favor small eccentricities (Kruckow et al. 2021). Such post–common-envelope systems most prominently accentuate the model differences discussed in Sect. 4. Namely, our fiducial model (S23b; Sect. 3) predicts that the orbit of accreting post–common-envelope systems will shrink; but only after their eccentricity has been excited to e ≈ 0.5. In contrast, the models based on Z21 and DD21 project that systems with low residual eccentricity after the plunge-in will circularize; and that the disk will widen their orbital separations instead of shrinking them. The latter, we note, is more consistent with the observed small eccentricities of post–common-envelope systems. Therefore, of the models considered, none provide a mechanism to reduce the orbital separation while maintaining consistency with the observed low eccentricities of post–common-envelope systems.

The previous considerations hinge on the assumption that the mass that can be accreted from a post–common-envelope circumbinary disk is sufficient to induce orbital change. The amount of residual envelope mass after the plunge-in phase, however, can vary greatly. For example, Passy et al. (2012) finds that a ≈ 1 M⊙ donor can retain more than 90% of the envelope. In contrast, Lau et al. (2022b) finds that the whole envelope of a ≈ 12 M⊙ donor can be ejected. In the case there is enough mass left to form a disk, Tuna & Metzger (2023) show that, for their chosen configuration, up to ΔM/M0 = 10% is accreted. This would be sufficient to pump (damp) binary eccentricity to the equilibrium band (circularity), but would not be adequate to significantly alter the orbital semi-major axis.

5.4. Post-AGB binaries

Circumbinary disks are also thought to form from the asymmetric outflows of AGB stars in binary systems (Mohamed & Podsiadlowski 2012; Saladino & Pols 2019), and are commonly observed around systems containing a post-AGB star (de Ruyter et al. 2006; van Winckel et al. 2009; Bujarrabal et al. 2013; Gallardo Cava et al. 2021; Kluska et al. 2022; Corporaal et al. 2023). Post-AGB systems are usually composed of a main-sequence star and a compact white dwarf or subdwarf companion, and have periods ranging from hundreds to thousands of days. A long-standing problem resides in the eccentricity of these systems, which for Porb < 2000 d should have been tidally circularized during the AGB phase (Pols et al. 2003; Izzard et al. 2010), however, they are found to have eccentricities up to e ≈ 0.6 (Jorissen et al. 1998; Oomen et al. 2018). Eccentricity pumping via a circumbinary disk has been invoked to explain this phenomenon (Dermine et al. 2013; Antoniadis 2014; Vos et al. 2015; Izzard & Jermyn 2023), although Rafikov (2016) argues against this possibility.

Moreover, Oomen et al. (2018) measure a bimodal distribution of eccentricities in post-AGB binaries with peaks around e = 0 and e = 0.25 (see also the discussion in Z21) which, if driven by interaction with a circumbinary disk, would suggest that near-circular accreting binaries tend to remain circular as in Z21 and DD21.

Our findings show that the mass accreted from the disk needs to be at least 5–10% the initial binary mass in order to drive significant orbit evolution, and the low masses (10−3 − 10−2 M⊙) contained in post-AGB disks (Bujarrabal et al. 2013) are likely insufficient. Moreover, the effects of disk winds and photoevaporation have not been included in the hydrodynamic models, but would likely amplify this issue.

5.5. Supermassive black hole binaries

Supermassive black hole binaries are expected to form as a result of galaxy mergers (Begelman et al. 1980). If gas is present in the nuclear region, it can settle into a rotationally supported circumbinary disk (Barnes 2002). The gas supply feeding the disk is self-regulated by the black hole accretion feedback, causing accretion to happen in short ≈105 yr bursts, during which the black holes can accrete mass close to the Eddington limit (Schawinski et al. 2015). Bursts are spaced out by quiescent phases, where the accretion luminosity can drop to 10−5 LEdd (Novak et al. 2011). Active and quiescent phases can alternate 100 to 10 000 times, comprising an Active Galactic Nucleus (AGN) phase lasting 107 − 109 yr (Martini 2004; Schawinski et al. 2015). Siwek et al. (2020) showed that, within a realistic cosmological framework, SMBH binaries with M < 108 M⊙ can double or triple their initial mass via disk accretion before they merge, whereas, as the initial mass of the binaries increases, their ability to significantly increase their mass diminishes. These results depend on the disk torque model employed (in this case based on the analytic work of Haiman et al. 2009) and would likely vary when using the prescriptions derived in this work. Nonetheless, disk accretion can significantly increase the mass of SMBH binaries, rendering them an ideal scenario for studying the long-term evolution of the orbit with our formalism.

However, the simulations underlying this work assume a disk aspect ratio of h/r = 0.1, which is notably higher than the typical values of h/r ≈ 10−2 − 10−3 expected for supermassive binaries accreting from stable, thin disks (see Appendix A); particularly relevant for LISA progenitors (105 − 107 M⊙). Moreover, observations of AGN variability point to a viscosity 0.01 ≤ α ≤ 0.03 (Starling et al. 2004; King et al. 2007), up to an order of magnitude lower than the value we assume. Simulating the low viscosity and/or low aspect ratio regime poses numerical challenges, but it may present with important morphological and evolutionary differences from the h/r = 0.1, α = 0.1 case (Tiede et al. 2020; Dittmann & Ryan 2022, 2024). Therefore, while our results are subject to change as one decreases viscosity and disk scale height, if a similar eccentricity attractor exists within this regime, we would expect a large population of supermassive binaries to maintain a relatively large eccentricity e ∼ 0.5, which could be manifest in electromagnetic searches. This eccentricity will decrease once gravitational-wave emission starts dominating the orbital evolution, but a detectable nonzero residual eccentricity should still be present when entering the LISA band (Cuadra et al. 2009; Roedig et al. 2012; Zrake et al. 2021; Garg et al. 2024). It is important to mention that we are only considering comparable-mass binaries (q > 0.1) in our discussion. The dynamics of SMBHs with much lower mass ratios are more similar to protoplanetary rather than stellar-binary systems (e.g., Armitage & Natarajan 2002).

Black hole binaries detectable through Pulsar Timing Arrays are more massive (above 108 − 109 M⊙) and have wider separations, and as such often lie outside the self-gravitating radius in standard α-disk-like solutions (Haiman et al. 2009). Thus, the models based on such solutions here presented (and Eq. (A.2)) are strictly speaking not directly applicable. However, alternative solutions can preserve the disk’s stability in such scenarios (Sirko & Goodman 2003; Gilbaum & Stone 2022), and gravitational stresses in an unstable disk may still be characterizable with an α-prescription (Rice et al. 2005). Moreover, Roedig et al. (2011) found that self-gravitating (but nonfragmenting) disks may display an equilibrium eccentricity around e ≈ 0.5, similar to the ones underpinning the models in this work. Such accretion around SMBH binaries in the PTA band could influence the power spectrum of the gravitational-wave background (Agazie et al. 2023a; EPTA Collaboration & InPTA Collaboration 2023; Reardon et al. 2023; Xu et al. 2023). For example, differential accretion can raise the median mass ratio of coalescing black holes and increase the power spectrum amplitude (Siwek et al. 2020). The power spectrum turnover at the lowest frequencies (which was tentatively observed Agazie et al. 2023b; EPTA Collaboration & InPTA Collaboration 2024) could also be accounted for with a population of eccentric binaries that shift gravitational wave emission to higher harmonic frequencies (although eccentricities higher than 0.5 are required to significantly alter the spectral shape) (Sesana et al. 2009; Sesana 2013; Chen et al. 2017) or with environment induced orbital decay (Kocsis & Sesana 2011; Agazie et al. 2023b); each of which would be consistent with sustained periods of circumbinary accretion.

6. Conclusion

In this study, we explore the long-term orbital evolution of binary systems interacting with circumbinary disks. Our model is based on the result of the hydrodynamic simulations by S23b that fully resolve the cavity region and the circumstellar disks, self-consistently accounting for gravitational and accretion forces. Our findings can be summarized as follows:

  • To alter significantly the orbital separation, the system needs to accrete a mass comparable to the total mass of the binary (Sects. 3.1 and 4.2).

  • If the initial eccentricity exceeds about 0.1, two consistent trends emerge (Sect. 3.1):

    • After accreting ≈10% of the initial mass, the eccentricity will stabilize at an equilibrium value e ≈ 0.5 that depends weakly on the mass ratio.

    • After accreting a mass comparable to the mass of the binary, the mass ratio will approach unity, and the orbital separation will shrink appreciably.

  • For initial eccentricities below approximately 0.1, our three models yield contrasting predictions. The model based on S23b suggests that the eccentricity will be excited to e ≈ 0.5, and then the orbit will shrink, while the models based on Z21 and DD21 predict that the orbit will circularize and then widen. The difference between the two predictions hinges on the different signs of the eccentricity derivative e′ when e < 0.1 (Sect. 4.1).

As a consequence of the above findings, a population of binaries that are interacting (or have recently interacted) with a thin, stable circumbinary disk will show some clear signatures in the orbital parameters distributions:

  • S23b: a correlation between mass ratio and eccentricity (Sect. 3.2).

  • Z21 and DD21: a gap in the eccentricity distribution at e ≈ 0.1 (Sect. 4.3).

We discuss the applicability of these results to various astrophysical scenarios, namely binary star formation (Sect. 5.1), interacting stars (Sect. 5.2), post–common-envelope systems (Sect. 5.3), post-AGB binaries (Sect. 5.4) and SMBH binaries (Sect. 5.5). We conclude that in almost all the considered scenarios, the systems are unable to accrete from the disk enough mass to significantly change the orbital separation. However, the interaction with the disk may still leave the above-mentioned signatures on the eccentricity distribution. A notable exception arises with SMBH binaries, as they can substantially increase mass through disk accretion. Consequently, the interaction with the disk will impact also the mass ratio and semi-major axis distribution, with implications for LISA and PTA observations. We note that the generality of these results is limited by our assumptions. We base our work on a grid of 2D simulations, while 3D effects are likely important (Duffell et al. 2024). We assume a fixed disk aspect ratio and viscosity prescription, while disks in nature can span a wide range of these parameters (see Fig. 9). Similarly, we ignore the effect of self-gravity, disk winds and disk-orbit inclination, which can all become important for specific applications. The vast parameter space that ensues from these variations has not been systematically explored and will require further investigation before these effects can be included in a formalism like the one we presented here.

Acknowledgments

R.V. thanks Norbert Langer, Shazrene Mohamed, Carlos Badenes, and all the participants of the MPA-Kavli Summer Program 2023 for the useful discussions and suggestions. A special thanks to Paul Ricker, for providing valuable insight on the (post-)common envelope phase, to Alice Somigliana, for the eye-opening discussion on protostellar disks, and to Deepika Bollimpalli, for answering far too many questions about AGN disks. RV and SdM are also grateful to Giuseppe Lodato and Alessia Franchini for discussions the numerical approach and 3D effects. The authors also thank the referee for a thorough and constructive report. D.D. and C.T. acknowledge support from the Danish Independent Research Fund through Sapere Aude Starting Grant No. 121587. J.C. acknowledges support from ANID, Millenium Nuclei NCN19_171 (NPF) and NCN2023_002 (TITANs). This research was supported in part by grant NSF PHY-1748958 to the Kavli Institute for Theoretical Physics (KITP). This project was part of the Kavli Summer Program in Astrophysics 2023, hosted at the Max Planck Institute for Astrophysics. We thank the Kavli Foundation and the MPA for their support. Software: NumPy (Harris et al. 2020), SciPy (Virtanen et al. 2020), Matplotlib (Hunter 2007).

References

  1. Adams, F. C., Ruden, S. P., & Shu, F. H. 1989, ApJ, 347, 959 [Google Scholar]
  2. Agazie, G., Anumarlapudi, A., Archibald, A. M., et al. 2023a, ApJ, 951, L8 [NASA ADS] [CrossRef] [Google Scholar]
  3. Agazie, G., Anumarlapudi, A., Archibald, A. M., et al. 2023b, ApJ, 952, L37 [NASA ADS] [CrossRef] [Google Scholar]
  4. Aly, H., Lodato, G., & Cazzoletti, P. 2018, MNRAS, 480, 4738 [NASA ADS] [Google Scholar]
  5. Andrew, S., Penoyre, Z., Belokurov, V., Evans, N. W., & Oh, S. 2022, MNRAS, 516, 3661 [NASA ADS] [CrossRef] [Google Scholar]
  6. Andrews, S. M., Chandler, C. J., Isella, A., et al. 2014, ApJ, 787, 148 [NASA ADS] [CrossRef] [Google Scholar]
  7. Antoniadis, J. 2014, ApJ, 797, L24 [NASA ADS] [CrossRef] [Google Scholar]
  8. Armitage, P. J., & Natarajan, P. 2002, ApJ, 567, L9 [NASA ADS] [CrossRef] [Google Scholar]
  9. Artymowicz, P. 1983, Acta Astron., 33, 223 [NASA ADS] [Google Scholar]
  10. Artymowicz, P., Clarke, C. J., Lubow, S. H., & Pringle, J. E. 1991, ApJ, 370, L35 [NASA ADS] [CrossRef] [Google Scholar]
  11. Backs, F., Poorta, J., Rab, Ch., et al. 2023, A&A, 671, A13 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  12. Barnes, J. E. 2002, MNRAS, 333, 481 [Google Scholar]
  13. Begelman, M. C., Blandford, R. D., & Rees, M. J. 1980, Nature, 287, 307 [Google Scholar]
  14. Belczynski, K., Taam, R. E., Kalogera, V., Rasio, F. A., & Bulik, T. 2007, ApJ, 662, 504 [Google Scholar]
  15. Blundell, K. M., Bowler, M. G., & Schmidtobreick, L. 2008, ApJ, 678, L47 [NASA ADS] [CrossRef] [Google Scholar]
  16. Bonnell, I. A., & Bate, M. R. 1994, MNRAS, 271, 999 [NASA ADS] [CrossRef] [Google Scholar]
  17. Bourne, M. A., Fiacconi, D., Sijacki, D., Piotrowska, J. M., & Koudmani, S. 2023, MNRAS, submitted [arXiv:2311.17144] [Google Scholar]
  18. Bujarrabal, V., Alcolea, J., Van Winckel, H., Santander-García, M., & Castro-Carrizo, A. 2013, A&A, 557, A104 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  19. Bujarrabal, V., Castro-Carrizo, A., Van Winckel, H., et al. 2018, A&A, 614, A58 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  20. Chen, S., Sesana, A., & Del Pozzo, W. 2017, MNRAS, 470, 1738 [CrossRef] [Google Scholar]
  21. Corporaal, A., Kluska, J., Van Winckel, H., et al. 2023, A&A, 674, A151 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  22. Cuadra, J., Armitage, P. J., Alexander, R. D., & Begelman, M. C. 2009, MNRAS, 393, 1423 [Google Scholar]
  23. de Ruyter, S., van Winckel, H., Maas, T., et al. 2006, A&A, 448, 641 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  24. de Vries, N., Portegies Zwart, S., & Figueira, J. 2014, MNRAS, 438, 1909 [NASA ADS] [CrossRef] [Google Scholar]
  25. Dempsey, A. M., Muñoz, D., & Lithwick, Y. 2020, ApJ, 892, L29 [NASA ADS] [CrossRef] [Google Scholar]
  26. Dermine, T., Izzard, R. G., Jorissen, A., & Van Winckel, H. 2013, A&A, 551, A50 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  27. Deroo, P., Acke, B., Verhoelst, T., et al. 2007, A&A, 474, L45 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  28. Di Stefano, R. 2020, MNRAS, 493, 1855 [CrossRef] [Google Scholar]
  29. Dittmann, A. J., & Ryan, G. 2021, ApJ, 921, 71 [NASA ADS] [CrossRef] [Google Scholar]
  30. Dittmann, A. J., & Ryan, G. 2022, MNRAS, 513, 6158 [NASA ADS] [Google Scholar]
  31. Dittmann, A. J., & Ryan, G. 2024, ApJ, 967, 12 [NASA ADS] [CrossRef] [Google Scholar]
  32. D’Orazio, D. J., & Duffell, P. C. 2021, ApJ, 914, L21 [CrossRef] [Google Scholar]
  33. D’Orazio, D. J., & Charisi, M. 2023, arXiv e-prints [arXiv:2310.16896] [Google Scholar]
  34. D’Orazio, D. J., Haiman, Z., & MacFadyen, A. 2013, MNRAS, 436, 2997 [Google Scholar]
  35. Dorozsmai, A., Toonen, S., Vigna-Gómez, A., de Mink, S. E., & Kummer, F. 2024, MNRAS, 527, 9782 [Google Scholar]
  36. Dubus, G., Taam, R. E., & Spruit, H. C. 2002, ApJ, 569, 395 [NASA ADS] [CrossRef] [Google Scholar]
  37. Duffell, P. C. 2016, ApJS, 226, 2 [NASA ADS] [CrossRef] [Google Scholar]
  38. Duffell, P. C., D’Orazio, D., Derdzinski, A., et al. 2020, ApJ, 901, 25 [NASA ADS] [CrossRef] [Google Scholar]
  39. Duffell, P. C., Dittmann, A. J., D’Orazio, D. J., et al. 2024, ApJ, submitted [arXiv:2402.13039] [Google Scholar]
  40. Dunham, M. M., Stutz, A. M., Allen, L. E., et al. 2014, The Evolution of Protostars: Insights from Ten Years of Infrared Surveys with Spitzer and Herschel (University of Arizona Press), 195 [Google Scholar]
  41. Dunhill, A. C., Cuadra, J., & Dougados, C. 2015, MNRAS, 448, 3545 [NASA ADS] [CrossRef] [Google Scholar]
  42. Dutrey, A., Guilloteau, S., & Simon, M. 1994, A&A, 286, 149 [NASA ADS] [Google Scholar]
  43. Dutrey, A., Guilloteau, S., Prato, L., et al. 1998, A&A, 338, L63 [NASA ADS] [Google Scholar]
  44. Dutrey, A., Di Folco, E., Beck, T., & Guilloteau, S. 2016, A&ARv, 24, 5 [Google Scholar]
  45. El-Badry, K., Rix, H.-W., Tian, H., Duchêne, G., & Moe, M. 2019, MNRAS, 489, 5822 [NASA ADS] [CrossRef] [Google Scholar]
  46. Elbakyan, V. G., Nayakshin, S., Vorobyov, E. I., Garatti, Caratti O. A., & Eislöffel, J., 2021, A&A, 651, L3 [Google Scholar]
  47. EPTA Collaboration & InPTA Collaboration (Antoniadis, J., et al.) 2023, A&A, 678, A50 [CrossRef] [EDP Sciences] [Google Scholar]
  48. EPTA Collaboration & InPTA Collaboration (Antoniadis, J., et al.) 2024, A&A, 685, A94 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  49. Farris, B. D., Duffell, P., MacFadyen, A. I., & Haiman, Z. 2014, ApJ, 783, 134 [NASA ADS] [CrossRef] [Google Scholar]
  50. Fiorellino, E., Tychoniec, Ł., Cruz-Sáenz de Miera, F., et al. 2023, ApJ, 944, 135 [NASA ADS] [CrossRef] [Google Scholar]
  51. Fleming, D. P., & Quinn, T. R. 2017, MNRAS, 464, 3343 [Google Scholar]
  52. Franchini, A., Sesana, A., & Dotti, M. 2021, MNRAS, 507, 1458 [NASA ADS] [Google Scholar]
  53. Franchini, A., Lupi, A., Sesana, A., & Haiman, Z. 2023, MNRAS, 522, 1569 [NASA ADS] [CrossRef] [Google Scholar]
  54. Fritsch, F. N., & Butland, J. 1984, SIAM J. Sci. Stat. Comput., 5, 300 [Google Scholar]
  55. Gagnier, D., & Pejcha, O. 2023, A&A, 674, A121 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  56. Gagnier, D., & Pejcha, O. 2024, A&A, 683, A4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  57. Gaia Collaboration (Arenou, F., et al.) 2023, A&A, 674, A34 [CrossRef] [EDP Sciences] [Google Scholar]
  58. Gallardo Cava, I., Gómez-Garrido, M., Bujarrabal, V., et al. 2021, A&A, 648, A93 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  59. Garg, M., Tiwari, S., Derdzinski, A., et al. 2024, MNRAS, 528, 4176 [Google Scholar]
  60. Gilbaum, S., & Stone, N. C. 2022, ApJ, 928, 191 [NASA ADS] [CrossRef] [Google Scholar]
  61. Gingold, R. A., & Monaghan, J. J. 1977, MNRAS, 181, 375 [NASA ADS] [CrossRef] [Google Scholar]
  62. Grant, S. L., Stapper, L. M., Hogerheijde, M. R., et al. 2023, AJ, 166, 147 [NASA ADS] [CrossRef] [Google Scholar]
  63. Guilloteau, S., Dutrey, A., & Simon, M. 1999, A&A, 348, 570 [NASA ADS] [Google Scholar]
  64. Gupta, A., Miotello, A., Manara, C. F., et al. 2023, A&A, 670, L8 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  65. Haiman, Z., Kocsis, B., & Menou, K. 2009, ApJ, 700, 1952 [CrossRef] [Google Scholar]
  66. Hamers, A. S., Glanz, H., & Neunteufel, P. 2022, ApJS, 259, 25 [NASA ADS] [CrossRef] [Google Scholar]
  67. Harris, C. R., Millman, K. J., Van Der Walt, S. J., et al. 2020, Nature, 585, 357 [NASA ADS] [CrossRef] [Google Scholar]
  68. Hartmann, L., & Kenyon, S. J. 1996, ARA&A, 34, 207 [NASA ADS] [CrossRef] [Google Scholar]
  69. Hunter, J. D. 2007, Comput. Sci. Eng., 9, 90 [NASA ADS] [CrossRef] [Google Scholar]
  70. Hwang, H.-C., Ting, Y.-S., & Zakamska, N. L. 2022, MNRAS, 512, 3383 [NASA ADS] [CrossRef] [Google Scholar]
  71. Izzard, R. G., & Jermyn, A. S. 2023, MNRAS, 521, 35 [NASA ADS] [CrossRef] [Google Scholar]
  72. Izzard, R. G., Dermine, T., & Church, R. P. 2010, A&A, 523, A10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  73. Jacobsen, S. K., Jørgensen, J. K., van der Wiel, M. H. D., et al. 2018, A&A, 612, A72 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  74. Jørgensen, J. K., van Dishoeck, E. F., Visser, R., et al. 2009, A&A, 507, 861 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  75. Jorissen, A., Van Eck, S., Mayor, M., & Udry, S. 1998, A&A, 332, 877 [NASA ADS] [Google Scholar]
  76. Kashi, A., & Soker, N. 2011, MNRAS, 417, 1466 [NASA ADS] [CrossRef] [Google Scholar]
  77. Kido, M., Takakuwa, S., Saigo, K., et al. 2023, ApJ, 953, 190 [NASA ADS] [CrossRef] [Google Scholar]
  78. King, A. R., Pringle, J. E., & Livio, M. 2007, MNRAS, 376, 1740 [NASA ADS] [CrossRef] [Google Scholar]
  79. Kley, W., & Nelson, R. P. 2012, ARA&A, 50, 211 [Google Scholar]
  80. Kluska, J., Van Winckel, H., Coppée, Q., et al. 2022, A&A, 658, A36 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  81. Kocsis, B., & Sesana, A. 2011, MNRAS, 411, 1467 [NASA ADS] [CrossRef] [Google Scholar]
  82. Kratter, K., & Lodato, G. 2016, ARA&A, 54, 271 [NASA ADS] [CrossRef] [Google Scholar]
  83. Kruckow, M. U., Neunteufel, P. G., Di Stefano, R., Gao, Y., & Kobayashi, C. 2021, ApJ, 920, 86 [NASA ADS] [CrossRef] [Google Scholar]
  84. Kuffmeier, M., Jensen, S. S., & Haugbølle, T. 2023, Eur. Phys. J. Plus, 138, 272 [NASA ADS] [CrossRef] [Google Scholar]
  85. Lacour, S., Biller, B., Cheetham, A., et al. 2016, A&A, 590, A90 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  86. Lau, M. Y. M., Hirai, R., González-Bolívar, M., et al. 2022a, MNRAS, 512, 5462 [NASA ADS] [CrossRef] [Google Scholar]
  87. Lau, M. Y. M., Hirai, R., Price, D. J., & Mandel, I. 2022b, MNRAS, 516, 4669 [NASA ADS] [CrossRef] [Google Scholar]
  88. Leigh, N. W. C., Toonen, S., Portegies Zwart, S. F., & Perna, R. 2020, MNRAS, 496, 1819 [NASA ADS] [CrossRef] [Google Scholar]
  89. Lesur, G. 2021, J. Plasma Phys., 87, 205870101 [CrossRef] [Google Scholar]
  90. Long, F., Andrews, S. M., Vega, J., et al. 2021, ApJ, 915, 131 [NASA ADS] [CrossRef] [Google Scholar]
  91. Lubow, S. H., & Artymowicz, P. 1996, Evolutionary Processes in Binary Stars. ASIC, 477, 53 [NASA ADS] [Google Scholar]
  92. Lucy, L. B. 1977, AJ, 82, 1013 [NASA ADS] [CrossRef] [Google Scholar]
  93. MacFadyen, A. I., & Milosavljević, M. 2008, ApJ, 672, 83 [Google Scholar]
  94. Manara, C. F., Natta, A., Rosotti, G. P., et al. 2020, A&A, 639, A58 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  95. Martin, R. G., Lepp, S., Zhang, B., Nixon, C. J., & Childs, A. C. 2024, MNRAS, 528, L161 [Google Scholar]
  96. Martini, P. 2004, Carnegie Observatories Astrophysics Series, 1, 169 [Google Scholar]
  97. Mathieu, R. D., Stassun, K., Basri, G., et al. 1997, AJ, 113, 1841 [NASA ADS] [CrossRef] [Google Scholar]
  98. Mignon-Risse, R., Oliva, A., González, M., Kuiper, R., & Commerçon, B. 2023, A&A, 672, A88 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  99. Miranda, R., Muñoz, D. J., & Lai, D. 2017, MNRAS, 466, 1170 [NASA ADS] [CrossRef] [Google Scholar]
  100. Mohamed, S., & Podsiadlowski, P. 2012, Balt. Astron., 21, 88 [Google Scholar]
  101. Moody, M. S. L., Shi, J.-M., & Stone, J. M. 2019, ApJ, 875, 66 [NASA ADS] [CrossRef] [Google Scholar]
  102. Moreno, M. M., Schneider, F. R. N., Röpke, F. K., et al. 2022, A&A, 667, A72 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  103. Muñoz, D. J., & Lithwick, Y. 2020, ApJ, 905, 106 [CrossRef] [Google Scholar]
  104. Muñoz, D. J., Springel, V., Marcus, R., Vogelsberger, M., & Hernquist, L. 2013, MNRAS, 428, 254 [CrossRef] [Google Scholar]
  105. Muñoz, D. J., Miranda, R., & Lai, D. 2019, ApJ, 871, 84 [CrossRef] [Google Scholar]
  106. Novak, G. S., Ostriker, J. P., & Ciotti, L. 2011, ApJ, 737, 26 [Google Scholar]
  107. Offner, S. S. R., Moe, M., Kratter, K. M., et al. 2023, ASP Conf. Ser., 534, 275 [NASA ADS] [Google Scholar]
  108. Oomen, G.-M., Van Winckel, H., Pols, O., et al. 2018, A&A, 620, A85 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  109. Pakmor, R., Springel, V., Bauer, A., et al. 2016, MNRAS, 455, 1134 [Google Scholar]
  110. Passy, J.-C., De Marco, O., Fryer, C. L., et al. 2012, ApJ, 744, 52 [NASA ADS] [CrossRef] [Google Scholar]
  111. Petzold, L. 1983, SIAM J. Sci. Stat. Comput., 4, 136 [CrossRef] [Google Scholar]
  112. Pols, O. R., Karakas, A. I., Lattanzio, J. C., & Tout, C. A. 2003, ASP Conf. Ser., 303, 290 [NASA ADS] [Google Scholar]
  113. Price, D. J., Cuello, N., Pinte, C., et al. 2018, MNRAS, 477, 1270 [Google Scholar]
  114. Pringle, J. E. 1991, MNRAS, 248, 754 [NASA ADS] [Google Scholar]
  115. Rafikov, R. R. 2016, ApJ, 830, 8 [NASA ADS] [CrossRef] [Google Scholar]
  116. Ragusa, E., Lodato, G., & Price, D. J. 2016, MNRAS, 460, 1243 [CrossRef] [Google Scholar]
  117. Ragusa, E., Alexander, R., Calcino, J., Hirsh, K., & Price, D. J. 2020, MNRAS, 499, 3362 [NASA ADS] [CrossRef] [Google Scholar]
  118. Ramírez-Tannus, M. C., Kaper, L., de Koter, A., et al. 2017, A&A, 604, A78 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  119. Ramírez-Tannus, M. C., Poorta, J., Bik, A., et al. 2020, A&A, 633, A155 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  120. Ramírez-Tannus, M. C., Backs, F., de Koter, A., et al. 2021, A&A, 645, L10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  121. Reardon, D. J., Zic, A., Shannon, R. M., et al. 2023, ApJ, 951, L6 [NASA ADS] [CrossRef] [Google Scholar]
  122. Reichardt, T. A., De Marco, O., Iaconi, R., Tout, C. A., & Price, D. J. 2019, MNRAS, 484, 631 [NASA ADS] [CrossRef] [Google Scholar]
  123. Riaz, R., Schleicher, D. R. G., Vanaverbeke, S., & Klessen, R. S. 2021, MNRAS, 507, 6061 [CrossRef] [Google Scholar]
  124. Rice, W. K. M., Lodato, G., & Armitage, P. J. 2005, MNRAS, 364, L56 [NASA ADS] [CrossRef] [Google Scholar]
  125. Ricker, P. M., & Taam, R. E. 2012, ApJ, 746, 74 [Google Scholar]
  126. Ricker, P. M., Timmes, F. X., Taam, R. E., & Webbink, R. F. 2019, IAU Symp., 346, 449 [Google Scholar]
  127. Roedig, C., Dotti, M., Sesana, A., Cuadra, J., & Colpi, M. 2011, MNRAS, 415, 3033 [NASA ADS] [CrossRef] [Google Scholar]
  128. Roedig, C., Sesana, A., Dotti, M., et al. 2012, A&A, 545, A127 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  129. Rosotti, G. P. 2023, New Astron. Rev., 96, 101674 [NASA ADS] [CrossRef] [Google Scholar]
  130. Röpke, F. K., & De Marco, O. 2023, Liv. Rev. Comput. Astrophys., 9, 2 [CrossRef] [Google Scholar]
  131. Saladino, M. I., & Pols, O. R. 2019, A&A, 629, A103 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  132. Sana, H., Ramírez-Tannus, M. C., de Koter, A., et al. 2017, A&A, 599, L9 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  133. Scaife, A. M. M. 2013, MNRAS, 435, 1139 [CrossRef] [Google Scholar]
  134. Schawinski, K., Koss, M., Berney, S., & Sartori, L. F. 2015, MNRAS, 451, 2517 [Google Scholar]
  135. Sesana, A. 2013, Class. Quant. Grav., 30, 224014 [NASA ADS] [CrossRef] [Google Scholar]
  136. Sesana, A., Vecchio, A., & Volonteri, M. 2009, MNRAS, 394, 2255 [NASA ADS] [CrossRef] [Google Scholar]
  137. Shakura, N. I., & Sunyaev, R. A. 1973, A&A, 24, 337 [NASA ADS] [Google Scholar]
  138. Shi, J.-M., Krolik, J. H., Lubow, S. H., & Hawley, J. F. 2012, ApJ, 749, 118 [Google Scholar]
  139. Shu, F. H., Tremaine, S., Adams, F. C., & Ruden, S. P. 1990, ApJ, 358, 495 [Google Scholar]
  140. Sirko, E., & Goodman, J. 2003, MNRAS, 341, 501 [NASA ADS] [CrossRef] [Google Scholar]
  141. Siwek, M. S., Kelley, L. Z., & Hernquist, L. 2020, MNRAS, 498, 537 [NASA ADS] [CrossRef] [Google Scholar]
  142. Siwek, M., Weinberger, R., Muñoz, D. J., & Hernquist, L. 2023a, MNRAS, 518, 5059 [Google Scholar]
  143. Siwek, M., Weinberger, R., & Hernquist, L. 2023b, MNRAS, 522, 2707 [NASA ADS] [CrossRef] [Google Scholar]
  144. Smallwood, J. L., Lubow, S. H., & Martin, R. G. 2022, MNRAS, 514, 1249 [NASA ADS] [CrossRef] [Google Scholar]
  145. Somigliana, A., Toci, C., Lodato, G., Rosotti, G., & Manara, C. F. 2020, MNRAS, 492, 1120 [NASA ADS] [CrossRef] [Google Scholar]
  146. Springel, V. 2010, MNRAS, 401, 791 [Google Scholar]
  147. Spruit, H. C., & Taam, R. E. 2001, ApJ, 548, 900 [NASA ADS] [CrossRef] [Google Scholar]
  148. Starling, R. L. C., Siemiginowska, A., Uttley, P., & Soria, R. 2004, MNRAS, 347, 67 [NASA ADS] [CrossRef] [Google Scholar]
  149. Stevenson, S., Vigna-Gómez, A., Mandel, I., et al. 2017, Nat. Commun., 8, 14906 [Google Scholar]
  150. Taam, R. E., & Spruit, H. C. 2001, ApJ, 561, 329 [CrossRef] [Google Scholar]
  151. Takakuwa, S., Saito, M., Saigo, K., et al. 2014, ApJ, 796, 1 [NASA ADS] [CrossRef] [Google Scholar]
  152. Tauris, T. M., & van den Heuvel, E. P. J. 2014, ApJ, 781, L13 [NASA ADS] [CrossRef] [Google Scholar]
  153. Testi, L., Natta, A., Manara, C. F., et al. 2022, A&A, 663, A98 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  154. Tiede, C., & D’Orazio, D. J. 2024, MNRAS, 527, 6021 [Google Scholar]
  155. Tiede, C., Zrake, J., MacFadyen, A., & Haiman, Z. 2020, ApJ, 900, 43 [NASA ADS] [CrossRef] [Google Scholar]
  156. Tiede, C., D’Orazio, D. J., Zwick, L., & Duffell, P. C. 2024, ApJ, 964, 46 [NASA ADS] [CrossRef] [Google Scholar]
  157. Tobin, J. J., Sheehan, P. D., Reynolds, N., et al. 2020, ApJ, 905, 162 [NASA ADS] [CrossRef] [Google Scholar]
  158. Tofflemire, B. M., Mathieu, R. D., Herczeg, G. J., Akeson, R. L., & Ciardi, D. R. 2017, ApJ, 842, L12 [CrossRef] [Google Scholar]
  159. Tokovinin, A. 2020, MNRAS, 496, 987 [NASA ADS] [CrossRef] [Google Scholar]
  160. Tokovinin, A., & Moe, M. 2020, MNRAS, 491, 5158 [NASA ADS] [CrossRef] [Google Scholar]
  161. Trani, A. A., Rieder, S., Tanikawa, A., et al. 2022, Phys. Rev. D, 106, 043014 [NASA ADS] [CrossRef] [Google Scholar]
  162. Tuna, S., & Metzger, B. D. 2023, ApJ, 955, 125 [NASA ADS] [CrossRef] [Google Scholar]
  163. Turpin, G. A., & Nelson, R. P. 2024, MNRAS, 528, 7256 [NASA ADS] [CrossRef] [Google Scholar]
  164. van Son, L. A. C., de Mink, S. E., Renzo, M., et al. 2022, ApJ, 940, 184 [NASA ADS] [CrossRef] [Google Scholar]
  165. van Winckel, H., Lloyd Evans, T., Briquet, M., et al. 2009, A&A, 505, 1221 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  166. Vigna-Gómez, A., MacLeod, M., Neijssel, C. J., et al. 2020, PASA, 37, e038 [Google Scholar]
  167. Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nat. Methods, 17, 261 [Google Scholar]
  168. Vos, J., Østensen, R. H., Marchant, P., & Van Winckel, H. 2015, A&A, 579, A49 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  169. Wang, H.-Y., Bai, X.-N., & Lai, D. 2023a, ApJ, 943, 175 [NASA ADS] [CrossRef] [Google Scholar]
  170. Wang, H.-Y., Bai, X.-N., Lai, D., & Lin, D. N. C. 2023b, MNRAS, 526, 3570 [NASA ADS] [CrossRef] [Google Scholar]
  171. Ward, W. R. 1997, ApJ, 482, L211 [NASA ADS] [CrossRef] [Google Scholar]
  172. Wei, D., Schneider, F. R. N., Podsiadlowski, P., et al. 2024, A&A, 688, A87 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  173. Weiser, A., & Zarantonello, S. E. 1988, Math. Comp., 50, 189 [CrossRef] [Google Scholar]
  174. Winter, A. J., & Haworth, T. J. 2022, Eur. Phys. J. Plus, 137, 1132 [NASA ADS] [CrossRef] [Google Scholar]
  175. Xu, H., Chen, S., Guo, Y., et al. 2023, Res. Astron. Astrophys., 23, 075024 [CrossRef] [Google Scholar]
  176. Young, M. D., & Clarke, C. J. 2015, MNRAS, 452, 3085 [NASA ADS] [CrossRef] [Google Scholar]
  177. Zrake, J., & MacFadyen, A. I. 2012, ApJ, 744, 32 [NASA ADS] [CrossRef] [Google Scholar]
  178. Zrake, J., Tiede, C., MacFadyen, A., & Haiman, Z. 2021, ApJ, 909, L13 [NASA ADS] [CrossRef] [Google Scholar]

Appendix A: Estimate of the disk aspect ratio

We can estimate the disk aspect ratio h/r of a thin accretion disk via equation 2.19 of Shakura & Sunyaev (1973), which, for parameters typical of an accreting stellar binary gives

h / r ≈ 0.12 ( M M ⊙ ) − 3 / 8 ( r AU ) 1 / 8 ( α 0.1 ) − 1 / 10 ( M ˙ 10 − 6 M ⊙ yr − 1 ) 3 / 20 , $$ \begin{aligned} \begin{aligned} h/r \approx \; 0.12&\left(\frac{M}{\mathrm{M}_{\odot }}\right)^{-3/8} \left(\frac{r}{\mathrm{AU}}\right)^{1/8} \left(\frac{\alpha }{0.1}\right)^{-1/10} \left(\frac{\dot{M}}{10^{-6} \,\mathrm{M}_{\odot }\,\mathrm{yr}^{-1}}\right)^{3/20}, \end{aligned} \end{aligned} $$(A.1)

where M is the mass of the central object (in this context, the total mass of the binary), r is the radial coordinate of the disk, to be taken as comparable to the size of the circumbinary cavity, α is the Shakura-Sunyaev viscosity parameter and Ṁ is the mass accretion rate.

When rewriting Equation A.1 for parameters typical of SMBHs, a significantly lower aspect ratio is obtained.

h / r ≈ 6.2 × 10 − 3 ( M 10 6 M ⊙ ) − 1 / 10 ( r 10 4 r g ) 1 / 8 ( α 0.1 ) − 1 / 10 ( L L Edd ) 3 / 20 ( η 0.06 ) − 3 / 20 , $$ \begin{aligned} \begin{aligned} h/r \approx \; 6.2 \times 10^{-3}&\left(\frac{M}{10^6 \,\mathrm{M}_{\odot }}\right)^{-1/10} \left(\frac{r}{10^4 \, r_g}\right)^{1/8} \\&\left(\frac{\alpha }{0.1}\right)^{-1/10} \left(\frac{L}{L_{\rm Edd}}\right)^{3/20} \left(\frac{\eta }{0.06}\right)^{-3/20}, \end{aligned} \end{aligned} $$(A.2)

where r g = 2 G M c 2 $ r_g=\frac{2GM}{c^2} $ is the Schwarzschild radius, L = η Ṁc2 is the disk luminosity, η is the gravitational energy release efficiency, and LEdd is the Eddington luminosity.

Appendix B: Additional figures

We describe two additional figures. Figure B.1 shows our interpolation of the results of DD21. As discussed in Sect. 2.3, the functions a′ and e′ by DD21 display ambiguous oscillations. They are likely caused by the precession of the disk cavity, which happens on a timescale of a few hundred orbits and therefore is not completely smoothed out by the ≈444 orbits window of the filter employed by DD21. We perform a discrete resampling of the curves, selecting the points in such a way as to smooth out the oscillations, while at the same time preserving the underlying trend and the points where the curves change sign. In the bottom panels, we show the residuals of the linear interpolation.

thumbnail Fig. B.1.

Interpolation of the results of DD21. In the top panel, a comparison between the values of a′ and e′ computed in DD21 (lighter solid lines) and the interpolating function used in this work (darker solid lines with dots). The bottom panels show the residuals.

Figure B.3 compares the results of the hydrodynamic simulations by S23b, Z21 and DD21 for q = 1. The general trend is similar among the three sets of simulations, the most relevant differences for the long-term evolution of the orbit are the discrepancy in the sign of e′ for e < 0.1 and the different location of the equilibrium eccentricity. The former difference implies that the predictions of the three models diverge for systems that start with e < 0.1 (Section 4.1). The latter determines a different speed in the semi-major axis evolution once systems have reached the equilibrium eccentricity (Section 4.2).

Figure B.2 shows how fast the evolution of the semi-major axis is for systems at the equilibrium eccentricity. Systems with small q < 0.2 widen their orbit, while for larger mass ratios the orbit always shrinks. The equilibrium eccentricity is reached after binaries accrete only about 10% of their initial mass (Section 3.1) and, after reaching it, they remain close to the equilibrium as the mass ratio steadily increases. Therefore, the mass ratio axis can be also interpreted as a time evolution axis. In this sense, a binary system that stars with q < 0.2 will initially increase its orbital separation, then it will transition to orbital shrinking as the mass ratio becomes larger than ≈0.3. The fastest orbital shrinking (a′≈ − 4) takes place when the binary passes from q ≈ 0.6. In the end, when the system reaches the attractor at q = 1, the rate of separation change settles to a′≈ − 0.7. We also show the corresponding values for Z21 and DD21 at q = 1.

thumbnail Fig. B.2.

Rate of separation change a′ computed at the equilibrium eccentricity, as a function of the mass ratio. We compare the models based on S23b (red), Z21 (green) and DD21 (blue). The blue cross indicates the value obtained by DD21 when considering small oscillations that a realistic system can have around the equilibrium eccentricity. Time moves from left to right as the binary accretes mass and the mass ratio approaches unity.

thumbnail Fig. B.3.

Results from the literature. Comparison of hydrodynamic simulations of an equal-mass binary and different eccentricities in a circumbinary disk from S23b (red), Z21 (green) and DD21 (blue). The dots correspond to the simulations and the solid lines correspond to our interpolation. The top panel shows the value of the derivative of the semi-major axis a′, while the bottom panel shows the derivative of the eccentricity e′. Vertical dotted lines mark the value of the equilibrium eccentricity, where e′ = 0.

Appendix C: Tables of derivatives

We report here for convenience the tables of derivatives used for the interpolation, as described in Sect. 2.

Table C.1.

Values of the semi-major axis derivative a′ from S23b.

Table C.2.

Values of the eccentricity derivative e′ from S23b.

Table C.3.

Values of the relative accretion rate λ from S23a.

Table C.4.

Values of the semi-major axis derivative a′ and eccentricity derivative e′ derived from the results of Z21 (columns 1-3) and DD21 (columns 4-7).

All Tables

Table 1.

Summary of the main parameters and results of the considered suites of hydrodynamic simulations.

In the text
Table C.1.

Values of the semi-major axis derivative a′ from S23b.

In the text
Table C.2.

Values of the eccentricity derivative e′ from S23b.

In the text
Table C.3.

Values of the relative accretion rate λ from S23a.

In the text
Table C.4.

Values of the semi-major axis derivative a′ and eccentricity derivative e′ derived from the results of Z21 (columns 1-3) and DD21 (columns 4-7).

In the text

All Figures

thumbnail Fig. 1.

Orbital evolution of binary systems under interaction with a circumbinary disk. The white lines follow the trajectories in the eccentricity–mass ratio plane. The thick black line traces the equilibrium eccentricity, that is the region where, for a fixed mass ratio q, the eccentricity derivative is e′ = 0. The background color shows our interpolated function approximating the derivative of the separation a′≡dlog a/dlog M. This quantity shows how fast the separation changes per unit of accreted mass. Regions where the orbit shrinks are shown in blue and regions where the orbit widens in red. The small black dots in the background correspond to the parameters used in the hydrodynamic simulations by S23b, on which the interpolation is based.
In the text
thumbnail Fig. 2.

Effect of accretion on the orbit of the binary. The top panel shows how much mass the systems need to accrete to reach the vicinity of the attractor point (white circle). The bottom panel shows the amount of accreted mass needed to reach the equilibrium band only (white region). The black line traces the equilibrium eccentricity. See also Fig. 1.
In the text
thumbnail Fig. 3.

Velocity at which systems travel across the eccentricity–mass ratio plane. The streamlines are the same as in Fig. 1, but here the color of the streamlines encodes how fast the mass ratio (top) and eccentricity (bottom) change per unit of accreted mass, quantified by the derivatives |q′| and |e′|. The black line traces the equilibrium eccentricity.
In the text
thumbnail Fig. 4.

Evolution of eccentricity (e) and semi-major axis (a) as a function of accreted mass. For each panel, the systems start from a different initial mass ratio (q0) that spans from 0.1 (top left) to 1.0 (bottom right). Within each panel, the initial eccentricity (e0) spans from 0.01 to 0.8. The trajectories are colored based on the mass change ΔM/M0, where M0 is the initial mass of the binary and ΔM is the amount of accreted mass. All the systems start along the vertical line a = a0, highlighted with a thick gray line.
In the text
thumbnail Fig. 5.

Evolution of the distribution of the orbital parameters for a population of binaries. Distribution of eccentricity (top), semi-major axis (middle), and mass ratio (bottom). Each color corresponds to a different amount of accreted mass, increasing from dark purple to light yellow.
In the text
thumbnail Fig. 6.

Population synthesis of binaries interacting with circumbinary disks in the e-q plane. Each panel corresponds to a different amount of accreted mass, increasing from top-left to bottom-right. The color scale indicates the number of systems falling into each bin. The gray band marks the position of the equilibrium eccentricity.
In the text
thumbnail Fig. 7.

Model comparison of the orbital evolution of equal-mass binaries. We show eccentricity (top) and semi-major axis (bottom) as a function of accreted mass. Each line represents a different system, colored according to the initial eccentricity, in the range 0.01 ≤ e0 ≤ 0.8. The black dashed lines in the bottom panels mark the analytical solutions obtained in Sect. 4.2 considering a system evolving at the equilibrium eccentricity.
In the text
thumbnail Fig. 8.

Model comparison of the evolution of the distributions of the orbital parameters for a population of equal-mass binaries. Distribution of eccentricity (top), semi-major axis (bottom), colored based on the amount of accreted mass. Each column shows one of the three models.
In the text
thumbnail Fig. 9.

Summary of the properties of circumbinary disks in astrophysical scenarios. Binary mass refers to the initial mass M0 of the central binary. The disk is fed if there is a mass inflow that replenishes the disk. The disk is self-gravitating if the self-gravity of the disk has a nonnegligible effect on the dynamics and structure of the disk. Disk mass Md is the instantaneous mass of the disk. Accretion rate is the mass accreted by the binary per unit time. Disk lifetime is the duration of a disk-mediated phase in the binary evolution, terminating when the disk is dispersed or the binary merges. Accreted mass is the amount of mass that can be accreted from the disk on to the binary within the disk lifetime. When the disk is not fed from external mass inflow, the accreted mass can be estimated as Md/M0, which assumes that all the mass in the disk is eventually accreted, and neglects outflows or other mechanisms of mass loss. h/r is the disk aspect ratio in the vicinity of the binary. α is the viscosity parameter of the disk. [1]Jørgensen et al. (2009), [2]Tobin et al. (2020), [3]Jacobsen et al. (2018), [4]Kido et al. (2023), [5]Fiorellino et al. (2023), [6]Takakuwa et al. (2014), [7]Manara et al. (2020), [8]Guilloteau et al. (1999), [9]Scaife (2013), [10]Dutrey et al. (2016), [11]Long et al. (2021), [12]Testi et al. (2022), [13]Grant et al. (2023), [14]Tuna & Metzger (2023), [15]Vos et al. (2015), [16]Izzard & Jermyn (2023), [17]Bujarrabal et al. (2018), [18]Kluska et al. (2022), [19]Bujarrabal et al. (2013), [20]Franchini et al. (2021), [21]Novak et al. (2011), [22]Schawinski et al. (2015), [23]Siwek et al. (2020), [24] estimated from Eq. (A.1) assuming r = 100 AU and α = 0.1, [25] estimated from Eq. (A.2) assuming r = 0.01 pc, α = 0.1 and η = 0.06, [26]Rosotti (2023), [27]Starling et al. (2004).
In the text
thumbnail Fig. B.1.

Interpolation of the results of DD21. In the top panel, a comparison between the values of a′ and e′ computed in DD21 (lighter solid lines) and the interpolating function used in this work (darker solid lines with dots). The bottom panels show the residuals.
In the text
thumbnail Fig. B.2.

Rate of separation change a′ computed at the equilibrium eccentricity, as a function of the mass ratio. We compare the models based on S23b (red), Z21 (green) and DD21 (blue). The blue cross indicates the value obtained by DD21 when considering small oscillations that a realistic system can have around the equilibrium eccentricity. Time moves from left to right as the binary accretes mass and the mass ratio approaches unity.
In the text
thumbnail Fig. B.3.

Results from the literature. Comparison of hydrodynamic simulations of an equal-mass binary and different eccentricities in a circumbinary disk from S23b (red), Z21 (green) and DD21 (blue). The dots correspond to the simulations and the solid lines correspond to our interpolation. The top panel shows the value of the derivative of the semi-major axis a′, while the bottom panel shows the derivative of the eccentricity e′. Vertical dotted lines mark the value of the equilibrium eccentricity, where e′ = 0.
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.