Extension of the Differential Algebraic Framework to the Three-Body Problem: A Constructive Approach to Unified Analytic Solutions in Celestial Mechanics

This paper establishes a constructive differential algebraic framework for the three-body problem by extending the nonlinear partial differential algebraic closure KNLPDE recently developed by the authors. We define the nonlinear ordinary differential algebraic closure K (3-body) NLODE as a differentially closed field extension constructed through a recursive adjunction process that incorporates: (i) fundamental solutions of linearized variational equations around two-body reference orbits, (ii) multi-index radical extensions Φ1/p, (iii) roots of unity ωp, and (iv) a constructive set of nonlinear special functions including Laplace coefficients, elliptic integrals, and hypergeometric functions arising from Fourier–Lagrange expansions. Within this closure, we prove that every analytic solution of the Newtonian three-body problem admits a unified local representation Mri (t) = r (0) i (t) +where r (0) im=1Φm (c, t) 1/pm ωkmpm ψm (t), i =1, 2, 3, is a two-body reference solution, ψm are nonlinear basis functions obtained via homotopy continuation and recursive Newton corrections from Keplerian variational modes, Φm are differential polynomials whose coefficients are the combinatorial correction coefficients Γ (2, 3) m, k, pm = 2 (Newtonian dominance), ωpm are primitive roots of unity, and km ∈ 0, 1 are branch selection indices uniquely determined by the initial data. We provide the first complete derivation of the combinatorial coefficients Γ (2, 3) m, k starting from the multivariate Faà di Bruno formula and the projection of quadratic nonlinearities onto the Fourier–Laplace basis of the linearized problem. Explicit formulas are given in terms of Beta functions, Laplace coefficients, and symmetry factors, and are shown to satisfy all necessary consistency conditions (normalization, causality, scaling, asymptotic exponential decay). Using these coefficients we: • Recover the classical Lagrange equilateral and Euler collinear solutions as special cases where Φm ≡ 0 or Φm reduces to an algebraic number of degree ≤ 5; • Derive the second-order amplitude–frequency relation for Lyapunov orbits around collinear libration points in the planar restricted three-body problem, reproducing and extending known expansions; • Obtain an algebraic proof of the Gascheau–Routh linear stability criterion for the Lagrange solution; • Construct the first explicit post-Newtonian (1PN) corrections to the combinatorial coefficients and use them to derive the existence condition for 1PN Euler solutions; • Propose a novel algebraic integrability criterion based on polynomial relations among the Φm polynomials, and prove that the three-body problem is nonintegrable for generic masses unless such a relation holds. All constructions are algorithmic: we present detailed pseudocode for precomputing Γ (2, 3) m, k (Algorithm 1) and for adaptive homotopy continuation with certified precision (Algorithm 2). Complexity analysis shows O (∥M∥ · ∥Bdmax ∥ + SM 3) operations, where M is the number of basis modes, S the homotopy steps, and Bdmax the set of nonlinear multi-indices. The method is validated on several benchmark problems (circular restricted three-body problem, Pythagorean three-body problem, equal-mass periodic orbit) using interval arithmetic; residuals below 10−12 are achieved with M = 16 and 128-bit precision. We further outline generalizations to the N -body problem (N ≥ 4), where the combinatorial coefficients become sums over unordered pairs and acquire a graphical interpretation, and to post-Newtonian approximations up to 3PN order. Finally, we discuss potential quantum speedup via amplitude estimation and HHL algorithm, and prove an algebraic criterion that characterizes integrable subcases solely through the structure of the Γ coefficients. This work demonstrates that the three-body problem, long considered unsolvable in closed form, admits explicit, constructive, and certifiable analytic solutions within a suitably enlarged differential algebraic universe. It unifies and transcends classical perturbation theories, provides new explicit results, and opens a systematic algebraic approach to nonlinear celestial dynamics.