Circumbinary Planet Stability Zone

Binary Star

Mass ratio μ = m₂/(m₁+m₂) 0.30

Binary eccentricity e_bin 0.20

Circumbinary Planet

Initial semi-major axis a_p (× a_bin) 3.00

Stability Prediction

Launch a planet to begin

Critical radius a_crit

– a_bin

Planet distance r_p

– a_bin

Max distance reached

– a_bin

Binary orbits elapsed

0.0

How it works

This is the restricted three-body problem: two stars of masses m₁, m₂ orbit their common barycenter on a fixed Keplerian ellipse (semi-major axis a_bin = 1, eccentricity e_bin), while a massless planet feels their combined, time-varying gravity and is integrated with real 4th-order Runge–Kutta:

a_planet = -G m1 (p-r1)/|p-r1|³ - G m2 (p-r2)/|p-r2|³
r1 = -μ R(t),  r2 = (1-μ) R(t)   [R(t) from Kepler's equation]

Whether a circumbinary ("P-type") orbit survives depends on how far out it starts. Holman & Wiegert (1999) fit the critical stability boundary from thousands of numerical integrations:

a_crit/a_bin = 1.60 + 5.10 e - 2.22 e² + 4.12 μ
             - 4.27 e μ - 5.09 μ² + 4.61 e² μ²

(valid for μ ∈ [0.1, 0.9], e ∈ [0, 0.8] — the same relation used to check the stability of real circumbinary planets like Kepler-16b.)

μ, e_bin sliders — set the binary's mass ratio and eccentricity; a_crit updates live from the formula above.
a_p slider + Launch — place the planet on a circular orbit at that radius and start the real N-body integration.
Snap to a_crit — sets a_p exactly to the predicted boundary so you can watch marginal cases slowly diverge.
Below a_crit, close encounters with the stars pump the planet's eccentricity chaotically until it is ejected; well above a_crit, the orbit stays quasi-periodic and bounded for as long as you run it.

Star 1 Star 2 Planet Critical radius a_crit

⚠ Couldn't load the 3D engineThree.js failed to load from the CDN. Check your connection and reload.